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Abstract 

The  Blast  Load  Simulator  (BLS)  at  the  U.S.  Army  Engineer  Research  and 
Development  Center  (ERDC)  has  been  utilized  for  applying  simulated 
blast  loads  on  structures  for  more  than  a  decade;  however,  the 
experimental  efforts  to  date  have  not  included  a  computational  modeling 
component.  The  Defense  Threat  Reduction  Agency  tasked  ERDC  to 
conduct  an  integrated  experimental/computational  program  to  evaluate 
several  first  principles  codes  (FPCs)  for  modeling  airblast  environments 
typical  of  those  encountered  in  the  BLS.  The  FPCs  considered  were  CTH, 
DYSMAS,  Loci/BLAST,  RAGE,  and  SHAMRC.  These  codes  afforded  a 
range  of  numerical  solvers  and  equation  of  state  options  for  the  modeling 
effort.  The  FPCs  were  evaluated  against  data  generated  in  the  BLS.  Three 
experimental  configurations  were  considered  for  the  validation.  These 
included  two  empty  configurations,  with  comparisons  drawn  against 
reflected  pressures  measured  at  the  downstream  end  wall.  The  third 
included  a  box-like  structure  in  the  flow  field.  The  structure  was 
instrumented  on  the  exposed  surfaces  to  provide  a  more  challenging 
scenario  for  the  model  comparisons.  The  purpose  of  the  integrated 
program  was  to  assess  computational  accuracy  and  cost  of  the 
aforementioned  FPCs  and  identify  any  shortcomings  in  the  physics 
modeling  and  areas  for  future  improvement.  Findings  from  the  integrated 
program  are  discussed  in  this  report. 
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1  Introduction 

1.1  Objective 

The  Blast  Load  Simulator  (BLS)  at  the  U.S.  Army  Engineer  Research  and 
Development  Center  (ERDC)  in  Vicksburg,  MS  has  been  utilized  for 
applying  simulated  blast  loads  on  structures  for  more  than  a  decade; 
however,  the  experimental  efforts  to  date  have  not  included  a 
computational  modeling  component.  The  Defense  Threat  Reduction 
Agency  (DTRA)  J9NTE  tasked  ERDC  to  conduct  an  integrated 
experimental  and  computational  program  to  evaluate  several  first 
principles  codes  (FPCs)  for  modeling  airblast  environments  typical  of 
those  encountered  in  the  BLS.  This  integrated  program  will  assess  both 
computational  accuracy  and  cost,  as  well  as  identify  any  shortcomings  in 
the  physics  modeling  and  areas  for  future  improvement. 

1.2  Purpose  and  scope 

DTRA  has  tasked  ERDC  to  conduct  an  integrated  experimental  and 
computational  program  to  evaluate  a  number  of  EPCs  commonly  used 
within  the  Department  of  Defense  (DoD)  community  for  modeling  air- 
blast.  The  codes  include  the  CTH  code  developed  by  Sandia  National 
Laboratories  (McGlaun  et  al.  1990;  Crawford  et  al.  2013b);  Dynamic 
System  Mechanics  Advanced  Simulation  (DYSMAS)  developed  by  the 
U.S.  Navy,  Naval  Surface  Warfare  Center,  Indian  Head  (Harris  et  al. 

2014);  Radiation  Adaptive  Grid  Eulerian  (RAGE)  developed  by  Leidos, 
formerly  Science  Applications  International  Corporation  (SAIC)  (Gittings 
et  al.  2008);  and  the  Second-Order  Hydrodynamic  Automatic  Mesh 
Refinement  Code  (SHAMRC)  developed  by  Applied  Research  Associates 
(Crepeau  et  al.  2012).  An  additional  code,  Loci/BLAST  (Luke  1999; 
Thompson  et  al.  2012b),  was  also  considered  in  the  evaluation.  The 
Loci/BLAST  code  is  under  development  at  the  Mississippi  State  University 
(MSU),  Center  for  Advanced  Vehicular  Systems.  The  DoD  codes  utilize 
Eulerian  solvers  of  various  types,  all  relying  on  a  structured  mesh.  The 
Loci/BLAST  code  utilizes  an  unstructured  mesh,  providing  a  well-rounded 
evaluation  of  capabilities  for  modeling  airblast.  There  are  other  codes 
capable  of  modeling  airblast;  however,  the  ones  chosen  are  thought  to 
have  a  fairly  broad  user  base  with  ongoing  support  on  a  number  of  DoD 
High  Performance  Computing  (HPC)  platforms. 
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The  FPCs  were  evaluated  against  data  generated  in  the  BLS.  Specifically, 
three  experimental  configurations  were  considered  for  the  current  model 
validation  effort.  These  included  two  empty  configurations  of  the  BLS, 
with  comparisons  drawn  against  reflected  pressures  measured  at  the 
downstream  end  wall.  The  third  case  included  a  box-like  structure  placed 
in  the  flow  field.  The  structure  was  instrumented  on  the  exposed  surfaces 
to  provide  a  more  challenging  scenario  for  the  model  comparisons. 

The  FPCs  evaluated  here  are  in  common  use  within  the  DoD.  They  utilize 
varied  numerical  approaches  for  solving  the  conservation  equations.  Thus, 
the  model  assessment  extends  beyond  just  matching  experimental  data  to 
include  the  viability  of  a  numerical  solver  to  model  this  class  of  airblast 
problem.  The  longer  term  goal  is  to  simulate  complex  three-dimensional 
(3-D)  environments,  which  include  a  number  of  obstacles  in  the  flow  field. 
It  is  important  to  understand  the  capabilities  (both  strengths  and 
weaknesses)  of  each  FPC  in  a  simplified  setting  before  tackling  more 
complicated  problems.  Although  several  of  the  BLS  configurations  could 
be  approximated  as  axisymmetric,  the  emphasis  for  the  FPC  comparisons 
was  on  3-D  modeling  as  this  provides  a  better  indication  of  a  code’s 
performance  for  handling  realistic  applications. 

This  report  focuses  on  the  computational  effort.  A  summary  of  the  BLS 
experimental  program  is  provided  herein  to  provide  background  for  the 
modeling  effort.  The  summary  should  provide  adequate  detail  of  the  three 
BLS  test  configurations  considered  and  insight  into  repeatability  of  the 
data.  Further  details  on  the  BLS  and  the  overall  experimental  program  can 
be  found  in  Dallriva  et  al.  (2016).  The  main  body  of  this  report  deals  with 
the  model  comparisons  and  assessment  of  FPC  capability.  One  should 
recognize  that  each  code  (and  solver)  is  unique.  The  analysts  conducting 
the  calculations  were  given  the  flexibility  to  model  the  BLS  tests  in  their 
own  way.  Further,  they  were  encouraged  to  explore  equation  of  state 
(EOS)  options,  conduct  mesh  resolution  studies,  and  investigate  varied 
approaches  for  modeling  the  BLS  driver  section.  An  appendix  is  devoted  to 
each  FPC  considered  to  allow  for  a  detailed  discussion  of  the  model  setup, 
the  reasoning  behind  choices  made  in  the  setup,  and  lessons  learned 
during  the  course  of  the  modeling.  The  main  body  of  this  report  focuses  on 
the  bigger  picture,  looking  for  overall  trends  and  areas  for  improvement 
needed  in  the  physics  modeling. 
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2  BLS  Overview  and  Test  Configurations 

2.1  General  description 

The  BLS,  located  at  ERDC,  is  housed  in  a  tunnel  structure  for  noise 
suppression  and  containment  as  depicted  in  Figure  i.  The  BLS  is 
composed  of  a  series  of  modular  sections  that  allow  flexibility  in 
prescribing  the  downstream  airblast  environment.  The  basic  configuration 
considered  for  this  work  is  depicted  in  Figure  2.  This  particular  setup  is 
referred  to  as  the  Government  Services  Administration  (GSA) 
configuration.  It  is  composed  of  (i)  a  constant  diameter  driver  section  with 
a  striker  mechanism  and  diaphragm,  (2)  a  continuously  vented  cone  (CVC) 
section  and  conical  expansion  section,  denoted  as  GSA  Ci,  (3)  a  second 
non-vented  conical  expansion  section,  denoted  as  GSA  C2,  (4)  two 
constant  diameter  transition  (TR)  sections,  denoted  as  GSA  TRi  and  TR2, 
(5)  a  constant  diameter,  telescoping  section,  denoted  as  the  GSA-Cascade, 
and  (6)  the  target  vessel.  All  of  the  sections  depicted  in  Figure  2  are 
circular.  It  should  be  reiterated  that  the  BLS  is  modular,  and  the  basic 
configuration  shown  here  is  only  one  of  many  possibilities.  A  general  over¬ 
view  of  the  BLS  is  provided  here.  Details  on  the  specific  configurations 
tested  will  be  discussed  in  subsequent  sections. 


Figure  1.  BLS  residing  within  tunnei  structure. 


Figure  2.  Basic  GSA  configuration. 
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The  driver  section  is  depicted  in  Figure  3.  The  driver  diameter  is  fixed, 
having  an  inside  diameter  of  16  in.  The  length  can  be  varied  from  18  to 
66  in.  in  6-in.  increments.  This  allows  flexibility  in  controlling  the  loading 
pulse  duration.  During  operations,  a  diaphragm  is  placed  on  the  down¬ 
stream  end  to  enclose  the  pressure  vessel.  The  driver  is  then  pressurized  to 
the  desired  test  (burst)  condition.  A  striker  mechanism  is  placed  outside  of 
the  diaphragm  and  subsequently  used  to  rupture  it.  Following  rupture,  the 
pressurized  gas  vents  into  the  downstream  BLS  sections.  The  striker 
assembly  (without  the  diaphragm)  is  depicted  in  Figure  4.  The  actual 
striker  mechanism  resides  outside  of  the  pressure  vessel.  The  flange  sec¬ 
tion  of  the  assembly  is  bolted  to  the  driver  with  the  diaphragm  on  the 
inside  to  seal  the  pressure  vessel  during  test  operations.  The  BLS  can 
accommodate  mixed  gases  (e.g.,  helium  and  air);  however,  for  the  tests 
reported  herein,  the  mixture  was  limited  to  100  percent  air  to  reduce 
uncertainty  in  the  EOS  modeling  by  the  FPCs. 

Figure  3.  Driver  section. 


(a)  Driver. 


Figure  4.  Striker  assembiy,  striker,  and  fiange. 


(a)  Exterior  view.  (b)  Interior  view. 
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The  rupture  of  the  diaphragm  results  in  the  generation  of  fragments  that 
travel  downstream.  These  fragments  can  adversely  impact  test  structures 
placed  within  the  BLS  as  well  as  damage  other  interior  components  and 
instrumentation.  A  grill  is  typically  included  in  the  test  setup  to  trap  frag¬ 
ments  and  mitigate  any  damage  to  downstream  components.  The  location 
of  the  grill  in  the  BLS  is  shown  in  Figure  5.  A  photograph  of  the  grill  is  pro¬ 
vided  in  Figure  6.  The  outside  diameter  of  the  grill  is  17-5/8  in.  The  hori¬ 
zontal  and  vertical  bars  in  the  grill  are  1/8-in.  thick  with  a  depth  of  1  in. 
The  spacing  between  the  horizontal  and  vertical  bars  is  1  and  4  in., 
respectively. 


Figure  5.  Grill  placement  in  the  BLS. 


The  grill  is  welded  at  the  end  of  the  first  conical  piece  from  the  diaphragm 
rupture  point.  The  outside  ring  fits  around  the  outside  of  the  first  cylindrical 
piece  at  the  first  vent  location.  The  grill  axial  length  Is  2-15/16  inch. 


Figure  6.  Photograph  of  grill. 


ERDC/GSL  TR-16-31 


6 


The  CVC  cross  section  depicted  in  Figure  5  is  composed  of  three  cylindrical 
sections,  referred  to  as  CVC-i,  -2,  and  -3.  CVC-i  has  a  length  and  inner 
diameter  of  31-7/8  and  16-1/16  in.,  respectively.  CVC-2  has  a  length  and 
inner  diameter  of  22-1/2  and  24  in.,  respectively.  CVC-3  has  a  length  and 
inner  diameter  of  22-1/2  and  32  in.,  respectively.  The  shell  thickness  for 
each  section  is  1/2  in.  The  first  and  second  sections  overlap  by  2-15/16  in. 
when  the  grill  is  attached  to  CYC-i.  There  are  cylindrical  vents  at  each 
overlap  of  the  CVC  sections.  This  is  the  primary  means  for  the  gas  blow¬ 
down  in  the  BLS. 

The  CVC  connects  to  an  expanding  conical  section  having  a  20-deg  cone 
angle.  The  length  of  the  expansion  section  is  79-13/32  in.  The  inside  diam¬ 
eter  of  this  cone  at  the  upstream  and  downstream  sides  is  40  and  68  in., 
respectively.  The  CVC  and  attached  cone  comprise  the  GSA  Ci  section.  The 
shell  thickness  is  1/2  in. 

The  GSA  C2  is  the  next  downstream  section.  It  also  has  a  20-deg  cone 
angle.  The  length  of  this  section  is  58-1/8  in.  The  inside  diameter  of  the 
C2  cone  at  the  upstream  and  downstream  sides  is  68  and  88-1/2  in., 
respectively.  The  shell  thickness  is  1/2  in. 

A  series  of  constant  diameter  sections  reside  downstream  of  GSA  C2. 

These  are  the  two  transitions  sections  (TRi  and  TR2)  and  the  Cascade.  The 
length  and  inner  diameter  of  each  of  the  transition  sections  is  46  and 
88-1/2  in.,  respectively.  The  shell  thickness  for  the  transition  sections  is 
1/2-in.  The  Cascade  is  a  telescoping  section  that  allows  for  a  variable 
length.  The  upstream  Cascade  section  has  a  length  and  inner  diameter  of 
58  and  88-1/2  in.,  respectively.  The  downstream  Cascade  section  has  a 
length  and  inner  diameter  of  52  and  91  in.,  respectively.  The  shell  thick¬ 
ness  for  each  is  1/2  in.  There  is  a  slight  radial  gap  between  the  Cascade 
sections  to  facilitate  the  telescopic  movement. 

The  target  vessel  at  the  far  end  of  the  BLS  serves  as  a  reaction  mass.  It  has 
a  flat  wall  on  one  end  (the  left  side  of  the  vessel  in  Figure  2)  that  has 
reflected  pressure  gauges  embedded  within  it.  This  is  commonly  referred 
to  as  either  the  calibration  plate  or  target  wall,  depending  on  the  nature  of 
the  test.  This  report  will  make  use  of  the  latter  term,  as  it  is  more  general. 
Reflected  pressure  data  on  the  target  wall  were  captured  in  all  of  the  tests 
outlined  in  this  report.  The  BLS  is  modular,  allowing  for  additional 
sections  to  be  included  in  the  configuration.  There  is  no  restriction  that  the 
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added  sections  be  circular.  In  fact,  two  of  the  configurations  discussed 
herein  had  square  cross  sections.  The  only  stipulation  is  that  the  transition 
of  the  blast  front  between  sections  be  smooth  and  have  a  negligible  effect 
on  the  waveform. 

The  expansion  and  transition  sections  downstream  of  the  CVC  are 
designed  to  help  facilitate  the  formation  of  a  planar  shock  front.  The 
amplitude  and  duration  of  the  pulse  are  largely  defined  by  the  burst 
pressure  and  the  driver  length.  An  ideal  waveform  is  one  that  has  an  initial 
peak  followed  by  an  exponential  decay.  In  practice,  the  character  of  the 
pulse  can  differ  from  the  ideal.  When  the  computational  effort  was  started, 
the  BLS  was  undergoing  a  number  of  design  iterations,  with  the 
waveforms  differing  slightly  from  the  ideal.  This  is  irrelevant  for  the 
modeling  effort,  as  the  FPCs  should  be  able  to  capture  the  form  and 
character  of  the  measured  waveform.  For  the  tests  reported  herein,  the 
burst  pressure  was  1,298  psig  for  all  cases.  This  ensured  that  the  initial 
conditions  were  consistent  for  the  FPC  comparisons. 

As  mentioned  earlier,  three  different  BLS  configurations  were  considered 
for  the  FPC  comparisons.  In  all  cases,  the  BLS  was  operating  in  a  closed 
configuration,  with  the  only  venting  occurring  in  the  CVC  section.  In  other 
words,  there  were  no  longitudinal  gaps  between  sections.  In  practice,  there 
was  some  venting  at  the  seals  between  sections;  however,  this  is  thought  to 
be  minimal,  with  venting  and  blow-down  controlled  mainly  by  the 
openings  in  the  CVC. 

The  BLS  configurations  considered  for  the  FPC  comparisons  are  outlined 
in  the  ensuing  discussion.  These  are  referred  to  by  their  base  configuration 
name,  or  simply  as  the  first,  second,  or  third  case  modeled.  It  is  important 
to  note  that  the  configuration  name  is  a  generic  term  and  does  not  refer  to 
any  specific  test  series.  The  naming  convention  used  here  is  strictly  for 
convenience  in  the  reporting  of  the  FPC  results. 

2.2  GSA  configuration,  caiibration  Tests  30a  and  30b 

The  first  case  modeled  considered  the  GSA  configuration  depicted  in 
Figure  2.  This  configuration  was  utilized  in  the  early  calibration  testing, 
which  varied  the  driver  length,  burst  pressure,  and  gas  mixture.  Tests  30a 
and  30b  were  considered  for  the  analysis.  The  setups  for  these  two  tests 
were  identical,  providing  some  indication  of  test  repeatability.  For  these 
tests,  the  driver  length  was  66  in.,  and  the  burst  pressure  was  1,298  psig. 
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The  driver  gas  was  composed  of  lOO-percent  air.  The  ambient  temperature 
was  nominally  77°F.  The  overall  length  of  the  Cascade  section  was 
98-3/16  in.,  resulting  in  an  overlap  of  11-13/16  in.  between  the  upstream 
and  downstream  sub-sections.  The  dimensions  of  all  other  BLS  sections 
are  as  discussed  in  the  general  overview. 


Gauges  were  embedded  in  the  target  wall  to  measure  the  reflected  pres¬ 
sure.  The  instrumentation  layout  is  shown  in  Figure  7.  Only  the  radial 
gauge  distances  were  available  at  the  time  the  analyses  were  started.  This 
was  not  an  issue  since  the  measured  data  exhibited  little  sensitivity  to 
gauge  location.  Thus,  the  analysts  could  specify  the  gauge  location  in  terms 
of  radial  distance  from  the  centerline,  as  the  measured  data  were  largely 
independent  of  azimuth. 


Figure  7.  Instrumentation  layout  on  target  wall  for  GSA  configuration. 


2.3  Calibration  tests  with  8-ftx8-ft  C2SQ  configuration 

The  second  case  modeled  was  an  extension  of  the  first  (see  Figure  8).  Two 
sections  were  added,  a  circle-to-square  transition  (C2SQ)  and  a  square 
section  (SQi).  A  square  target  vessel  with  an  8-ft  by  8-ft  target  mated  with 
section  SQi.  Dimensioned  drawings  for  these  additions  are  provided  in 
Figures  9  through  11.  An  interior  view  looking  towards  the  target  vessel  is 
provided  in  Figure  12.  Section  C2SQ  did  not  mate  perfectly  with  the 
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Cascade  resulting  in  an  overlap  at  the  radial  extents.  The  largest  overlap 
on  the  left-hand  side  was  1-1/4  in-  (with  respect  to  the  view  towards  the 
target  vessel).  Similarly,  the  largest  overlap  on  the  right-hand  side  was 
2-3/8  in. 

A  side-on  pressure  gauge  PR  was  placed  in  the  center  of  GSA  TR2.  The 
axial  location  is  depicted  in  Figure  8.  The  gauge  resided  on  the  centerline 
of  the  BLS  and  was  placed  23  ft  8-7/16  in.  downstream  of  the  driver.  A 
typical  mounting  system  for  the  side-on  pressure  gauge  is  shown  in  Fig¬ 
ure  13.  In  addition,  the  target  wall  was  instrumented  to  measure  the 
reflected  pressure  across  the  surface.  The  instrumentation  layout  is  shown 
in  Figure  14. 


Figure  8.  BLS  8-ftx8-ft  C2SQ  configuration. 
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Figure  9.  Circie-to-square  transition  section  C2SQ. 
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Figure  10.  Square  section  SQl. 


Figure  11.  Square  transition  for  target  vessei  face. 
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Figure  12.  View  towards  target  vessei. 


Figure  13.  Typicai  mounting  system  for  the  side-on  pressure  gauge. 


(a)  Mount  to  sidewaii  (b)  Gauge,  side  view  (c)  Gauge,  end- 

on  view 
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Figure  14.  Instrumentation  iayouton  target  waii. 


2.4  BLS  8-ft><8-ft  C2SQ  configuration  with  box 

For  the  third  case  modeled,  the  BLS  configuration  was  identical  to  that  of 
the  second  with  the  exception  of  a  box  structure  placed  in  section  C2SQ 
(see  Figure  15).  An  interior  view  is  provided  in  Figure  16,  where  the  viewer 
is  looking  towards  the  target  vessel.  The  box  was  positioned  along  the 
centerline  of  the  BLS  (see  Figure  17).  The  length,  width,  and  height  of  the 
box  were  13, 13,  and  18  in.,  respectively.  It  sat  on  a  pedestal,  so  that  its 
base  was  level  with  the  Cascade  floor.  The  box  was  instrumented  with 
pressure  gauges  to  characterize  the  time-varying  flow  around  it.  The  box 
was  instrumented  on  all  exposed  surfaces,  as  depicted  in  Figures  18 
through  20.  In  the  nomenclature,  PBF  refers  to  the  front  side  gauges;  PBB 
refers  to  the  back  side;  PBT  refers  to  the  top;  PBL  refers  to  the  left  side; 
and  PBR  refers  to  the  right  side.  The  orientation  of  the  left  and  right  sides 
are  defined  looking  from  the  target  vessel  towards  the  driver,  i.e.,  opposite 
of  the  gunner’s  point  of  view.  There  was  reduced  instrumentation  on  the 
target  wall,  with  only  gauges  P3  and  P9  being  active  (see  Figure  21). 
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Figure  15.  BLS  8  ftx8  ft  C2SQ  configuration  with  box. 


Figure  17.  Box  iocation  and  overaii  dimensions. 
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Figure  18.  Box  instrumentation,  front  and  back  sides. 


(a)  Front  side 

Figure  19.  Box  instrumentation,  top. 


Figure  20.  Box  instrumentation,  ieft  and  right  sides. 


(a)  Left  side 


(b)  Right  side 
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Figure  21.  Active  gauges  on  target  waii. 
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3  Modeling  Approach  and  Comparisons 

3.1  Overview 

The  FPCs  chosen  for  the  BLS  modeling  effort  were  CTH,  DYSMAS,  Loci/ 
BLAST,  RAGE,  and  SHAMRC.  With  the  exception  of  Loci/BLAST,  the 
codes  listed  are  in  common  use  within  the  DoD  for  modeling  airblast 
applications.  The  DoD  codes  utilize  Eulerian  solvers  of  various  types,  all 
relying  on  a  structured  mesh.  The  Loci/BLAST  code  differs  in  that  it  uti¬ 
lizes  an  unstructured  mesh.  There  are  other  codes  capable  of  modeling 
airblast;  however,  the  ones  chosen  are  thought  to  have  a  fairly  broad  user 
base  with  ongoing  support  on  a  number  of  DoD  HPC  platforms. 

The  chosen  codes  represent  an  array  of  numerical  approaches  for  model¬ 
ing  airblast  applications.  Our  goal  was  to  first  assess  the  capability  of  each 
for  modeling  the  airblast  environment  encountered  within  the  BLS.  This 
entailed  a  close  look  at  the  strengths  and  weaknesses  of  each  code  for 
modeling  this  application,  and  determining  any  shortfalls  in  the 
underlying  physics  models  that  should  be  addressed.  The  assessment  also 
considered  the  computational  cost  and  feasibility  of  a  particular  FPC  for 
modeling  a  BLS  event.  Although  the  first  two  cases  can  be  approximated  as 
axisymmetric,  the  emphasis  was  on  g-D  modeling  as  this  provides  a  better 
indication  of  a  code’s  performance  for  handling  realistic  applications. 

Thus,  3-D  calculations  were  conducted  for  all  cases  discussed  in  this 
report;  however,  planes  of  symmetry  were  permitted  across  the  centerline 
axis  of  the  BLS.  The  use  of  symmetry  conditions  was  at  the  modeler’s 
discretion. 

The  results  from  each  FPC  were  compared  against  the  experimental  data 
to  assess  accuracy  and  trends.  These  comparisons  were  qualitative  in 
nature,  as  no  in-depth  uncertainty  analysis  was  conducted.  It  should  be 
noted  that  repeat  experiments  were  performed  for  all  of  the  cases  mod¬ 
eled.  In  general,  the  experimental  data  were  highly  repeatable.  Some 
differences  were  noted  later  in  time;  however,  this  should  not  be 
surprising  as  the  timescale  spanned  out  to  the  6o-  to  lOO-msec  range.  The 
comparisons  focused  on  the  time  of  arrival  (TOA),  pressure,  and  impulse. 
The  emphasis  was  not  just  on  capturing  the  peaks,  but  also  considered  the 
ability  of  the  codes  to  capture  the  time-dependent  nature  of  the  waveform. 
The  history  was  important,  as  it  provided  insight  into  the  complex  wave 
interaction  that  occurred  later  in  time. 
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All  of  the  calculations  discussed  in  this  report  were  open,  i.e.,  the  analysts 
had  access  to  the  experimental  data  prior  to  conducting  the  calculation. 
The  goal  of  the  computational  effort  was  to  have  an  honest  and  open 
assessment  of  each  code’s  capability.  Regular  teleconferences  were  con¬ 
ducted  to  foster  communication  between  the  analysts  and  encourage  open¬ 
ness.  Without  this  openness,  one  cannot  get  a  clear  assessment  of  the 
strengths  and  weaknesses  of  the  individual  codes,  nor  a  path  forward  for 
future  model  development.  Future  efforts  will  move  toward  blind  calcula¬ 
tions,  where  the  modeler  does  not  have  access  to  the  measured  data 
beforehand. 

Table  i  provides  a  brief  description  of  the  FPCs  used  in  the  BLS  modeling. 
Further  details  on  the  numerical  solver  and  features  unique  to  each  code 
are  outlined  in  the  appendices.  An  appendix  is  devoted  to  each  FPC  con¬ 
sidered  to  allow  for  a  detailed  discussion  of  the  model  setup,  the  reasoning 
behind  choices  made  in  the  setup,  and  lessons  learned  during  the  course  of 
the  modeling.  The  discussion  for  each  code  is  in  alphabetical  order  to 
avoid  any  perception  of  favoritism.  Similarly,  the  order  of  the  appendices 
follows  this  alphabetical  listing.  Thus,  CTH  is  discussed  in  Appendix  A, 
DYSMAS  in  Appendix  B,  etc. 

The  emphasis  in  this  section  is  on  the  bigger  picture,  which  involves  over¬ 
all  comparisons  of  FPC  performance.  Performance  entails  an  assessment 
of  accuracy  (i.e.,  ability  to  replicate  the  measured  data)  as  well  as  compu¬ 
tational  cost  in  conducting  3-D  analyses.  The  overarching  goal  is  to  assess 
the  strengths  and  weaknesses  of  each  of  the  FPCs  and  identify  potential 
areas  for  improvement  needed  in  the  physics  modeling. 

Before  proceeding  with  the  comparisons,  it  is  worthwhile  to  first  outline 
the  modeling  challenges  posed  by  the  BLS  setup.  These  challenges  arise 
due  to  the  disparate  length  and  time  scales  that  will  be  encountered  in  the 
modeling.  The  thin  shell  of  the  BLS  structure  (0.5-in.)  drives  much  of  the 
resolution  requirement.  In  general,  one  needs  several  cells  across  the 
thickness  of  this  thin  structure  to  adequately  capture  it  in  the  model. 
Meeting  this  resolution  requirement  can  lead  to  excessively  long  run  times 
-  due  to  the  resulting  large  number  of  cells  and  small  time-step  -  that  are 
intractable  for  3-D  analyses.  It  is  not  just  a  matter  of  resolving  the  shell 
thickness,  but  also  the  need  to  adequately  resolve  the  shell  along  the  full 
length  of  the  BLS  (about  40  ft).  Furthermore,  there  are  additional 
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Table  1.  Brief  description  of  FPCs  used  in  BLS  modeling. 


FPC 

Brief  Description 

CTH 

Structured  Eulerian  code  with  two-step  solution  approach  involving  a 
Lagrangian  step  followed  by  a  remap.  Explicit  finite  volume 
approximation  applied  to  conservation  equations  for  Lagrangian  step 
on  a  staggered  mesh.  Remap  based  on  second-order-accurate  van 

Leer  scheme.  Code  supports  wide  range  of  EOS  options  and 
sophisticated  material  strength  models. 

DYSMAS 

Coupled  Euler/Lagrange  code.  Structured  Eulerian  code  Gemini  solves 
the  Eulerian  problem  using  a  time-split,  second  order  Godunov  method 
based  on  a  Monotonic  Upstream  Centered  Scheme  for  the 
Conservations  Laws  (MUSCL).  Lagrangian  problem  solved  using  the 
ParaDyn  explicit  Lagrangian  finite  element  code.  DYSMAS  loosely 
couples  the  two  solution  domains. 

Loci/BLAST 

Unstructured,  fully  Eulerian,  CFD  code  that  uses  a  Harten-Lax-van 
Leer-Enfeldt  (HLLE)  approximate  Riemann  solver  to  compute  fluxes  at 
cell  faces  and  least-squares  reconstruction  coupled  with  a  nodal 
Baldwin-Barth  limiter  to  provide  second-order  spatial  accuracy.  Explicit 
time  integration  performed  using  a  two-stage  second-order  Total 
Variation  Diminishing  (TVD)  Runge-Kutta  scheme. 

RAGE 

Structured,  fully  Eulerian,  CFD  code  that  uses  a  two-shock  or  single¬ 
intermediate-state  approximate  Riemann  solver  with  an  alternating- 
direction-explicit  technique  in  the  shock  capturing.  Robust  AMR  that  is 
continuous  in  both  space  and  time  with  adaption  triggered  by  changes 
in  state  and/or  first-order  truncation  error. 

SHAMRC 

Structured  Eulerian  code  that  uses  a  two-step  approach  involving  a 
Lagrangian  step  followed  by  a  remap.  Conservation  equations  solved 
using  a  modified  Lax-Wendroff  finite  difference  approximation  on  a 
staggered  mesh.  Remap  based  on  second-order  accurate  van  Leer 
scheme  with  an  alternating  direction  integration  scheme. 

complicated  thin  structures  such  as  the  catch  grill  that  may  require 
increased  resolution  if  included  in  the  model. 

The  time  scale  is  also  an  issue.  In  general,  it  takes  20  to  25  msec  for  the 
shock  front  to  traverse  the  full  length  of  the  BLS  (varies  depending  on  the 
configuration).  The  time  duration  needed  to  fully  capture  the  loading  pulse 
is  on  the  order  of  50  msec.  Time-steps  are  likely  on  the  microsecond 
timescale.  Thus,  it  will  take  a  considerable  number  of  cycles  to  conduct  the 
calculation  to  completion.  So  one  should  expect  that  any  BLS  calculation 
will  be  highly  resolved  and  long  running. 
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3.2  GSA  configuration,  caiibration  Tests  30a  and  30b 

The  modeling  began  with  the  GSA  configuration,  which  represented  the 
simplest  of  the  three  cases.  The  geometry  can  easily  be  approximated  as 
axisymmetric  if  the  flange  cross-bars  and  grill  are  omitted  from  the  model. 
The  modelers  took  advantage  of  this  fact  by  conducting  two-dimensional 
cylindrical  (2-DC)  analyses  early  on.  This  allowed  them  to  investigate 
modeling  options  such  as  mesh  resolution,  EOS  selection,  hydrodynamic 
treatment  options,  etc.  in  a  2-DC  setting  that  allowed  high  turnaround 
calculations.  The  early  2-DC  modeling  served  as  a  foundation  for  all 
follow-on  3-D  analyses.  All  discussion  of  the  early  2-DC  modeling  is 
provided  in  the  appendices.  The  reporting  in  this  section  focuses  solely  on 
the  3-D  analysis  results. 

The  approach  taken  by  each  modeler  differed.  This  is  reflected  in  the 
summary  provided  in  Table  2.  The  table  summarizes  parameters  pertinent 
for  ensuing  code-to-code  comparisons.  These  include  the  symmetry 
conditions,  mesh  resolution,  EOS  selected  for  the  air,  approximate  total 
cell  count,  and  inclusion  of  the  grill  in  the  model.  Further  details  on  prob¬ 
lem  setup  and  assumptions  can  be  found  in  the  appendices.  In  most  cases, 
the  resolution  was  fixed  within  the  region  of  interest,  i.e.,  the  immediate 
region  encompassing  the  BLS;  however,  there  are  cases  where  the  resolu¬ 
tion  varied  within  this  domain.  Loci/BLAST  utilizes  an  unstructured  mesh. 
The  mesh  was  composed  of  tetrahedral  elements.  The  nominal  resolution 
within  the  air  space  was  2.5  cm,  with  a  finer  resolution  employed  near  the 
vents  and  grill.  The  RAGE  calculation  utilized  AMR.  The  resolution 


Table  2.  Summary  of  FPC  problem  setup  for  GSA  configuration. 


FPC 

Model 

Geometry 

Mesh 

Resolution 

(cm) 

EOS  for  Air 

Approximate 
No.  of  Cells 
(millions) 

Grill 

Included? 

CTH 

Half¬ 

symmetry 

0.5 

SESAME 

395 

No 

DYSMAS 

Fully  3-D 

1.0 

Ideal  Gas 

136 

Yes 

Loci/BLAST 

Fully  3-D 

1.3  -  2.5 

Ideal  Gas 

52 

Yes 

RAGE 

Quarter- 

symmetry 

0.5  (AMR) 

Ideal  Gas 

13  -  279 

No 

SHAMRC 

Half¬ 

symmetry 

0.5 

Doan-Nickel 

1024 

No 
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reported  is  the  size  of  the  finest  cell  in  the  mesh.  The  use  of  AMR  also 
resulted  in  a  continual  growth  of  the  problem  size.  The  cell  counts 
reported  in  the  table  are  for  the  beginning  and  end  of  the  calculation. 

Most  modelers  only  ran  a  single  calculation;  however,  for  SHAMRC,  two 
variations  were  considered  for  this  first  case.  The  first  was  a  coarser- 
meshed  i.o-cm  resolution  calculation  that  was  fully  3-D.  The  second  was  a 
medium-meshed  0.5-cm  resolution  calculation  that  employed  a  single 
plane  of  symmetry.  The  latter  was  chosen  for  the  subsequent  comparisons 
due  to  its  higher  resolution.  As  outlined  in  Appendix  E,  there  was  little 
difference  in  results  for  the  two  calculations. 

Comparisons  of  the  measured  and  calculated  pressure  histories  at  select 
gauges  located  on  the  target  wall  are  provided  in  Figure  22.  The  figure  only 
contains  comparisons  for  three  gauges,  P3,  P4,  and  PE.  The  measured 
(and  calculated)  results  were  generally  consistent  across  the  target  wall, 
and  one  need  only  consider  a  subset  to  assess  code  performance.  The 
radial  distance  from  the  BLS  centerline  for  these  gauges  was  2.25, 13.2, 
and  30  in.,  respectively.  Thus,  the  gauges  selected  span  the  full  range  of 
radial  offsets  across  the  target  wall.  Unless  stated  otherwise,  the  calculated 
results  were  time-shifted  so  that  their  arrival  times  coincided  with  the 
measured  arrival.  This  was  done  to  ease  comparisons  of  the  waveforms.  In 
all  cases,  the  calculated  TOA  lagged  the  measured  arrival  time.  The  time- 
shifts  applied  to  the  target  wall  gauges  for  each  FPC  are  outlined  in  Table 
3.  An  example  of  the  non-shifted  results  for  gauge  PE  is  provided  in 
Figure  23  to  illustrate  the  relative  difference  in  arrival  times  calculated  by 
each  of  the  FPCs.  It  can  be  difficult  to  interpret  results  for  an  individual 
code  from  the  overlaid  data  in  Figure  22.  Comparisons  on  an  individual 
basis  are  provided  for  gauge  PE  in  Figure  24.  These  comparisons  are  rep¬ 
resentative  of  those  noted  at  the  other  gauge  locations. 

In  general,  one  notes  excellent  overall  agreement  with  the  measured  data 
for  CTH,  DYSMAS,  and  Loci/BLAST.  The  general  character  of  the  calcu¬ 
lated  waveform  for  SHAMRC  is  excellent,  with  the  exception  of  the 
overprediction  of  the  initial  peak  pressure.  The  overprediction  is  on  the 
order  of  20  to  25  percent.  This  initial  overprediction  leads  to  a 
corresponding  overprediction  in  the  integrated  impulse.  RAGE  exhibits 
similar  issues,  but  to  a  much  greater  degree  with  further  overprediction  of 
the  initial  peak  pressure  (on  the  order  of  35  to  40  percent).  It  is  worth 
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noting  that  these  same  overpredictions  were  noted  in  the  2-DC  RAGE 
calculations 


Figure  22.  FPC  comparisons  with  measured  data,  GSA  configuration. 


(a)  Gauge  P3.  (b)  Gauge  P4. 


(c)  Gauge  PE. 

Tabie  3.  Time-shifts  appiied  to  GSA 
configuration  caicuiations. 


FPC 

Target  Wall  Gauges 
(msec) 

CTH 

-4.59 

DYSMAS 

-2.16 

Loci/BLAST 

-3.42 

RAGE 

-3.15 

SHAMRC 

-3.72 

Impulse  (psi-msec) 
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Figure  23.  Early  time  comparison 
for  gauge  PE  without  time-shifting. 


Time  (msec) 


Figure  24.  Individual  FPC  comparisons  for  gauge  PE. 


(a)  CTH. 


(b)  DYSMAS. 


Time  (msec) 


(c)  Loci/BLAST. 


(d)  RAGE. 


Impulse  (psi-msec) 


ERDC/GSL  TR-16-31 


23 


Figure  24.  (Continued). 


(see  Appendix  D).  Both  convergence  and  sensitivity  studies  were  con¬ 
ducted;  however,  they  did  not  shed  insight  as  to  the  source  of  the  severe 
overpredictions  noted  with  RAGE.  Due  to  time  constraints,  it  was  decided 
to  move  forward  with  the  3-D  RAGE  calculations  with  the  expectation  that 
it  would  overpredict  both  the  initial  peak  pressure  and  integrated  impulse. 

All  codes  exhibit  increased  noise  in  the  solution  at  gauge  P3.  This  gauge  is 
near  to  the  centerline  of  the  BLS  (2.25-in.  offset).  The  BLS  walls  and  any 
symmetry  planes  in  a  calculation  act  as  perfect  reflectors,  and  it  is  conjec¬ 
tured  that  this  idealization  induces  ringing  in  the  solution  due  to  the  axi- 
symmetric  nature  of  the  problem  setup.  The  ringing  is  still  evident  in  the 
fully  3-D  calculations  of  DYSMAS  and  Loci/BLAST,  but  to  a  slightly  lesser 
extent  as  compared  to  those  using  a  plane  of  symmetry.  The  ringing 
behavior  is  clearly  exaggerated  in  the  RAGE  calculation,  and  appears  to  be 
pervasive  throughout  the  calculation.  The  ringing  in  the  solution  at  the 
various  gauge  locations  dampens  out  with  increasing  radial  offset  from  the 
centerline. 

The  ringing  behavior  noted  at  gauge  P3  is  evident  in  the  measured  data, 
but  to  a  far  lesser  degree.  This  is  illustrated  in  Figure  25  where  the 
measured  data  for  the  three  gauges  from  Testgoa  are  overlaid.  In  reality, 
the  walls  of  the  BLS  are  not  perfect  reflectors  either  due  to  surface 
roughness  or  perturbations  in  the  shell  thickness  along  the  length.  These 
imperfections  offset  the  timing  of  reflections  and  effectively  cancel  out 
much  of  the  shock  superposition  that  is  observed  in  the  calculations. 
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Figure  25.  Gauge  comparisons 
for  Test  30a. 


Time  (msec) 


3.3  Calibration  tests  with  8-ftx8-ft  C2SQ  configuration 

The  8-ftx8-ft  C2SQ  configuration  is  an  extension  of  the  first  case.  It 
involves  the  addition  of  the  C2SQ  and  SQi  sections  to  the  standard  GSA 
configuration.  All  modeling  of  this  case  was  conducted  in  3-D.  Table  4 
summarizes  the  problem  setup  associated  with  each  FPC.  Further  details 
on  problem  setup  and  assumptions  can  be  found  in  the  appendices.  A 
comparison  of  the  measured  and  calculated  pressure  histories  at  select 
gauges  is  provided  in  Figure  26.  The  figure  contains  comparisons  for  three 
gauges  on  the  target  wall  (P6,  P5,  and  Pi)  along  with  the  side-on  pressure 
gauge  PR.  The  radial  distances  of  the  target  wall  gauges  from  the  BLS  cen¬ 
terline  were  10.625,  24.5,  and  46.72  in.  for  gauges  P6,  P5,  and  Pi,  respec¬ 
tively.  Thus,  the  gauges  selected  provide  a  good  indication  of  the  influence 
of  radial  offset  on  the  data.  The  measured  (and  calculated)  data  were  fairly 
consistent  across  the  target  wall,  indicating  that  the  incoming  shock  front 
was  relatively  planar.  All  calculated  results  were  time-shifted  to  match  the 
measured  arrival  time.  The  time  shifting  was  done  to  ease  comparisons  of 
the  waveform.  The  time-shifts  applied  for  each  FPC  calculation  are  out¬ 
lined  in  Table  5.  It  can  be  difficult  to  interpret  results  for  an  individual 
code  from  the  overlaid  data  in  Figure  26.  Comparisons  on  an  individual 
basis  are  provided  for  gauge  Pi  in  Figure  27.  These  comparisons  are  repre¬ 
sentative  of  those  noted  at  the  other  gauge  locations.  Similarly,  individual 
comparisons  for  the  side-on  pressure  gauge  PR  are  provided  in  Figure  28. 

In  general,  one  notes  excellent  overall  agreement  with  the  measured  data 
at  the  target  wall  for  CTH,  DYSMAS,  and  Loci/BLAST.  It  is  apparent  that 
both  the  initial  peak  and  rate  of  decay  are  captured  well  in  the  first  pulse. 
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Table  4.  Summary  of  FPC  problem  setup  for  8-ftx8-ft 
C2SQ  configuration. 


FPC 

Model 

Geometry 

Mesh 

Resolution 

(cm) 

EOS  for 

Air 

Approximate 
No.  of  Cells 
(millions) 

Grill 

Included? 

CTH 

Half-symmetry 

0.5 

SESAME 

473 

No 

DYSMAS 

Fully  3-D 

1.0 

Ideal  Gas 

166 

Yes 

Loci/BLAST 

Fully  3-D 

1.3  -  2.5 

Ideal  Gas 

53 

Yes 

RAGE 

Half-symmetry 

0.5  (AMR) 

Ideal  Gas 

39  -  509 

No 

SHAMRC 

Half-symmetry 

0.5 

Doan- 

Nickel 

1088 

Yes 

Figure  26.  FPC  comparisons  with  measured 
data,  8-ftx8-ft  C2SQ  configuration. 


Time  (msec) 

(a)  Gauge  PI. 


Time  (msec) 

(b)  Gauge  P5. 
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(c)  Gauge  P6. 


(d)  Gauge  PR. 


ERDC/GSL  TR-16-31 


26 


Table  5.  Time-shifts  applied  to  8-ftx8-ft 
C2SQ  configuration  calculations. 


FPC 

Target  Wall  Gauges 
(msec) 

Side-on  Gauge  PR 
(msec) 

CTH 

-4.80 

-4.80 

DYSMAS 

-1.15 

0 

Loci/BLAST 

-3.67 

-3.96 

RAGE 

-2.93 

-3.79 

SHAMRC 

-3.92 

-4.06 

SHAMRC  falls  into  the  next  tier  and  is  a  very  close  second.  The  general 
character  of  the  calculated  waveform  for  SHAMRC  is  excellent;  however, 
there  is  an  overprediction  of  the  peak  pressure.  The  overprediction  is 
transient  in  nature  but  sufficient  to  induce  an  earlier  rise  in  the  impulse, 
and  hence,  an  overprediction  of  the  final  value.  If  not  for  these  transients, 
SHAMRC  would  have  performed  equally  well.  Overall,  RAGE  results 
compared  poorly.  The  initial  peak  was  significantly  overpredicted,  which 
results  in  a  significant  overprediction  of  the  impulse  (on  the  order  of 
40  percent  at  the  time  the  calculation  was  terminated). 

The  comparisons  for  the  side-on  pressure  gauge  PR  are  mixed.  All  codes 
capture  the  form  and  character  of  the  initial  peak  associated  with  the 
pulse;  however,  a  number  of  them  have  drop-offs  during  the  initial  decay. 
This  is  clearly  evident  for  CTH,  DYSMAS,  and  RAGE.  All  of  the  codes 
capture  the  later  time  spike;  however,  the  timings  differ.  CTH,  Loci/ 
BLAST,  RAGE,  and  SHAMRC  all  capture  the  timing  and  amplitude  of  the 
first  spike  reasonably  well.  Both  SHAMRC  and  RAGE  appear  to  capture 
the  timing  of  successive  spikes  exceedingly  well,  albeit  there  are  differ¬ 
ences  in  amplitude  that  are  most  notable  for  SHAMRC.  Beginning  at  about 
20  msec,  there  is  significant  variability  in  the  impulse  derived  from  the 
measured  data.  This  is  likely  the  result  of  the  interface  (contact  surface) 
between  the  driver  and  drier  air  impinging  on  the  gauge.  This  contact 
surface  produces  turbulence  that  affects  the  degree  of  repeatability  among 
the  tests. 


Overpressure  (psi)  Overpressure  (psi) 
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Figure  27.  Individual  FPC  comparisons  for  gauge  PI. 


Time  (msec) 

(a)  CTH. 


Time  (msec) 

(c)  Loci/BLAST. 


(b)  DYSMAS. 


(d)  RAGE. 


& 

> 

O 


Time  (msec) 


(e)  SHAMRC. 


Impulse  (psi-msec)  Impulse  (psi-msec) 


Overpressure  (psi)  Overpressure  (psi) 


ERDC/GSL  TR-16-31 


28 


Figure  28.  Individual  FPC  comparisons  for  side-on  gauge  PR. 
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3.4  BLS  8-ftx8-ft  C2SQ  configuration  with  box 

The  8-ftx8-ft  C2SQ  with  box  configuration  is  an  extension  of  the  second 
case.  It  involves  the  addition  of  an  instrumented  box  in  the  C2SQ  section. 
All  modeling  of  this  case  was  conducted  in  3-D.  Table  6  summarizes  the 
problem  setup  associated  with  each  FPC.  Further  details  on  problem  setup 
and  assumptions  can  be  found  in  the  appendices.  All  calculated  results 
were  time-shifted  to  match  the  measured  arrival  time.  The  time-shifts 
applied  for  each  FPC  calculation  are  outlined  in  Table  7.  A  comparison  of 
the  measured  and  calculated  pressure  histories  at  select  gauges  on  the  box 
and  target  wall  are  provided  in  Figures  29  through  33.  The  comparisons 
are  limited  to  one  gauge  on  each  face  of  the  box,  as  well  as  a  single  gauge 
on  the  target  wall.  The  comparisons  shown  are  representative  of  those 
noted  at  the  other  gauge  locations  on  each  box  face  and  the  target  wall. 

Comparison  of  results  for  this  third  case  is  much  more  subjective  given  the 
increased  complexity  of  the  BLS  configuration.  When  viewing  the  data  for 
box  gauges  (measured  and  calculated),  one  notes  two  distinct  peaks.  The 
first  is  associated  with  the  original  incident  wave  propagating  in  the  down¬ 
stream  direction.  The  second  is  associated  with  the  reflected  wave  that  is 
traveling  upstream  after  reflecting  off  of  the  target  wall.  One  should  expect 
the  latter  to  be  the  more  challenging  of  the  two  to  capture  computa¬ 
tionally.  Overall,  one  notes  excellent  correlation  between  the  measured 
and  calculated  results  for  Loci/BLAST  and  SHAMRC.  Here,  the  complete 
waveform  is  replicated  exceedingly  well,  with  only  a  slight  overprediction 
for  the  initial  peak  pressure.  The  next  tier  is  composed  of  CTH  and 
DYSMAS.  Both  codes  capture  the  general  form  and  character  of  the  wave¬ 
form;  however,  disparities  arise  in  the  rate  of  decay  associated  with  each 
pulse.  This  is  more  of  an  issue  for  CTH,  where  there  is  a  sudden  drop-off  in 
the  pressure  history  following  arrival  of  the  reflected  shock.  This  results  in 
a  noticeable  underprediction  in  the  impulse  at  later  times.  Overall,  RAGE 
performed  poorly.  In  particular,  RAGE  tends  to  overpredict  the  peak 
pressure  for  the  reflected  shock,  with  a  much  slower  decay  than  observed 
experimentally.  The  net  effect  is  a  substantial  overprediction  in  impulse  at 
later  times.  With  the  exception  of  RAGE,  all  codes  provided  reasonable 
correlation  with  the  measured  data  for  the  target  wall  gauges. 
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Table  6.  Summary  of  FPC  problem  setup  for  8-ftx8-ft 
C2SQ  with  box  configuration. 


FPC 

Model 

Geometry 

Mesh 

Resolution 

(cm) 

EOS  for 

Air 

Approximate 
No.  of  Cells 
(millions) 

Grill 

Included? 

CTH 

Half-symmetry 

0.5 

SESAME 

473 

No 

DYSMAS 

Fully  3-D 

1.0 

Ideal  Gas 

166 

Yes 

Loci/BLAST 

Fully  3-D 

1.3  -  2.5 

Ideal  Gas 

69.9 

Yes 

RAGE 

Half-symmetry 

0.5  (AMR) 

Ideal  Gas 

39  -  478 

No 

SHAMRC 

Half-symmetry 

0.5 

Doan- 

Nickel 

1088 

Yes 

Table  7.  Time-shifts  applied  to  8-ftx8-ft  C2SQ  with 
box  configuration  calculations. 


FPC 

Box  Gauges 
(msec) 

Target  Wall  Gauge 
(msec) 

CTH 

-4.11 

-4.11 

DYSMAS 

-1.43 

-0.63 

Loci/BLAST 

-3.18 

-2.78 

RAGE 

-2.93 

-2.49 

SHAMRC 

-3.50 

-3.30 

Overpressure  (psi)  Overpressure  (psi)  Overpressure  (psi) 
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Figure  29.  FPC  comparisons  for  front  face  gauge  PBFl. 
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Figure  30.  FPC  comparisons  for  top  face  gauge  PBT2. 
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Figure  31.  FPC  comparisons  for  right  face  gauge  PBRl. 
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Figure  32.  FPC  comparisons  for  back  face  gauge  PBBl. 
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Figure  33.  FPC  comparisons  for  target  waii  gauge  P9. 
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3.5  Qualitative  assessment  of  accuracy 

A  qualitative  assessment  of  code  accuracy  will  be  discussed  here  and 
subsequently  used  to  rank  the  performance  of  the  FPCs  evaluated  in  this 
study.  The  accuracy  assessment  is  subjective  and  relies  on  visual 
comparisons  between  the  measured  and  calculated  results  outlined  in  the 
previous  sections.  Although  no  quantitative  metrics  were  employed,  the 
qualitative  assessment  is  quite  useful  as  it  provides  an  indication  of 
whether  or  not  a  particular  code  can  even  be  applied  to  a  BLS-type  of 
problem.  Furthermore,  the  three  BLS  cases  considered  allow  the 
assessment  to  cover  a  range  of  physical  conditions.  This  is  important  to 
gain  insight  into  whether  or  not  a  code  can  provide  reasonable  predictions 
of  airblast  behavior,  rather  than  a  “postdictive”  response. 

The  code  rankings  are  based  on  tiers.  Each  code  is  designated  to  fall  into 
one  of  three  tiers,  with  codes  falling  into  the  first  tier  being  the  best  per¬ 
formers.  It  was  decided  to  develop  separate  rankings  for  the  three  different 
BLS  configurations  that  were  modeled.  Further,  it  was  decided  to 
subdivide  rankings  for  the  last  two  BLS  configurations  into  interior  and 
target  wall  gauges.  For  the  8-ftx8-ft  C2SQ  configuration,  the  side-on 
pressure  gauge  PR  is  the  interior  gauge.  For  the  8-ftx8-ft  C2SQ 
configuration  with  a  box,  the  gauges  residing  on  the  box  are  all  lumped 
into  the  category  of  interior  gauges.  A  side-by-side  comparison  of  the 
overall  waveform  was  made  for  each  of  the  codes  in  the  development  of  the 
rankings.  It  is  not  possible  to  assign  an  ordered  ranking  for  each  code  (1  to 
5)  due  to  differences  noted  throughout  the  full  time-history  of  the 
waveform.  Some  codes  might  capture  the  initial  peak  well  but  not  the 
second  peak.  Conversely,  others  might  overpredict  the  initial  peak  but 
track  the  measured  waveform  almost  exactly  thereafter. 

Qualitative  rankings  for  code  accuracy  are  provided  in  Table  8.  These  are 
based  on  the  author’s  opinion,  and  the  reader  is  free  to  disagree.  RAGE 
was  found  to  perform  poorly  overall  and  was  placed  into  the  third  tier.  It  is 
conjectured  that  there  is  an  issue  with  the  ideal-gas  EOS  implementation 
that  should  be  addressed.  This  is  discussed  in  detail  in  Appendix  D.  It  was 
clear  from  the  early  2-DC  scoping  calculations  that  RAGE  would 
significantly  overpredict  the  initial  peak  pressure  and  subsequently  the 
peak  impulse.  Given  the  time  constraints  on  the  project,  it  was  decided  to 
move  forward  with  the  3-D  calculations  with  the  expectation  that  the 
comparisons  would  be  poor.  This  was  somewhat  of  a  less  than  ideal 
approach  for  conducting  the  analysis  with  RAGE. 
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Table  8.  Qualitative  rankings  for  code  accuracy  by  tiers. 


FPC 

GSA 

Configuration 
Target  Wall 

8-ftx8-ftC2SQ 

Configuration 

8-ftx8-ftC2SQ 
Configuration  with  Box 

Target 

Wall 

Side-on 
Gauge  PR 

Target 

Wall 

Box 

Gauges 

CTH 

1 

1 

2 

1 

2 

DYSMAS 

1 

1 

2 

2 

1 

Loci/BLAST 

1 

1 

1 

1 

1 

RAGE 

3 

3 

3 

3 

3 

SHAMRC 

2 

2 

1 

1 

1 

For  the  first  case  (GSA  configuration),  the  first  tier  codes  were  CTH, 
DYSMAS,  and  Loci/BLAST.  CTH  seems  to  capture  the  overall  waveform 
the  best  of  the  three,  with  DSYMAS  and  Loci/BLAST  generally  under¬ 
predicting  the  arrival  time  of  the  second  pulse.  SHAMRC  was  placed  in  the 
second  tier  due  to  the  notable  overprediction  of  the  initial  peak  pressure. 
Otherwise,  it  replicated  the  overall  waveform  exceedingly  well.  It  should 
probably  be  considered  a  close  second. 

For  the  second  case,  a  similar  ranking  is  noted  for  the  target  wall  gauges 
for  all  of  the  same  reasons  as  discussed  for  the  first  case.  Again,  SHAMRC 
is  a  close  second,  only  because  of  the  notable  overprediction  in  the  initial 
peak  pressure.  Otherwise,  the  calculated  waveform  correlates  well  with  the 
measured  data.  Greater  disparities  between  the  calculated  and  measured 
data  are  noted  for  the  side-on  gauge  PR.  Both  Loci/BLAST  and  SHAMRC 
capture  the  general  structure  of  the  measured  waveform  as  well  as  the 
timing  of  the  later  time  first  spike.  SHAMRC  actually  captures  the  timing 
of  the  follow-on  spikes  exceedingly  well,  albeit  the  amplitudes  of  the  spikes 
are  significantly  over-predicted.  Both  CTH  and  DYSMAS  were  placed  in 
the  second  tier  due  to  the  drop-off  in  pressure  following  the  initial  peak.  It 
is  possible  that  the  tracers  resided  in  mixed  cells,  which  could  have 
affected  the  results.  This  might  explain  the  errant,  but  isolated,  differences 
between  the  calculated  and  measured  waveforms. 

There  was  a  slight  switch  for  the  target  wall  ranking  in  the  third  case, 
where  CTH,  Loci/BLAST,  and  SHAMRC  fall  into  the  first  tier  and 
DYSMAS  the  second.  The  slight  overprediction  of  the  initial  peak  pressure 
and  underprediction  of  the  arrival  time  for  the  second  pulse  resulted  in 
DYSMAS  being  ranked  in  the  second  tier.  One  should  note  that  neither  the 
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CTH  nor  Loci/BLAST  calculations  were  conducted  long  enough  to  capture 
the  second  pulse.  It  is  possible  that  they  may  have  also  predicted  an  earlier 
arrival  for  the  second  pulse.  Here,  DYSMAS  should  be  considered  a  close 
second. 

DSYMAS,  Loci/BLAST,  and  SHAMRC  all  fall  into  the  first  tier  for  compar¬ 
isons  against  the  box  gauges.  Both  Loci/BLAST  and  SHAMRC  replicate  the 
measured  waveforms  exceedingly  well,  at  least  for  the  duration  of  the  cal¬ 
culation.  There  is  slightly  more  disparity  between  the  DYSMAS-calculated 
waveform  and  that  measured,  which  is  generally  related  to  the  timing  of 
the  third  and  subsequent  pulses.  For  CTH,  there  are  noticeable  drop-offs 
in  the  calculated  pressure  following  the  arrival  of  the  reflected  shock 
(second  pulse).  The  drop-offs  appear  to  be  pervasive  for  the  front,  top,  and 
right  face  gauges. 

3.6  Computational  cost 

As  mentioned  earlier,  the  assessment  of  code  performance  entails  both 
accuracy  and  computational  cost.  The  discussion  so  far  has  focused  on 
accuracy,  albeit  the  assessment  was  subjective  and  relied  on  visual 
comparisons  between  the  measured  and  calculated  results.  The  focus  here 
is  on  assessing  the  cost  of  conducting  a  calculation  for  a  BLS-type  of 
problem.  The  end  goal  is  to  move  beyond  the  BLS  and  address  modeling 
large-scale,  real-world  environments  composed  of  multiple  structures  in 
the  blast  field.  Thus,  an  assessment  of  computational  cost  associated  with 
the  BLS  modeling  should  provide  an  indicator  if  an  FPC  can  move  forward 
for  production  computing  involving  more  realistic  applications. 

A  cost  comparison  for  each  case  is  provided  in  Tables  9  through  11.  The 
tables  contain  the  nominal  mesh  resolution  for  reference,  the  approximate 
number  of  cells,  the  end  time  of  the  calculation,  the  number  of  processors 
utilized,  and  the  wall-clock  hours  (or  central  processing  unit  (CPU)  time  if 
reported).  The  total  number  of  processor  hours  is  shown  in  the  last  col¬ 
umn.  This  metric  is  the  summation  of  the  number  of  processors  times  the 
wall-clock  time  (or  CPU  if  specified)  over  all  restarts.  The  summation  is 
important  as  it  captures  changes  in  the  resource  requirements.  This  is  best 
exemplified  by  a  RAGE  calculation,  where  the  number  of  processors  uti¬ 
lized  grew  as  the  calculation  progressed  forward  in  time.  CPU  times  are 
preferred,  since  they  exclude  the  expense  of  reading  and  writing  data; 
however,  the  wall-clock  times  were  typically  reported  by  the  modelers. 
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One  should  also  note  that  the  FPC  calculations  were  conducted  on  differ¬ 
ent  HPC  platforms.  Processor  speeds  and  memory  capacity  varied  for  each 
of  the  platforms.  Thus,  there  can  be  no  “apples-to-apples”  comparisons  as 
the  speed  of  computation  on  these  platforms  differ.  Regardless,  it  is 
possible  to  perform  a  cost  assessment  on  an  order-of-magnitude  basis  as 
all  of  the  HPC  architectures  utilized  the  current  state-of-the-art 
technology.  The  reader  should  be  able  to  recognize  which  FPCs  provide  a 
reasonable  turn-around  time  for  an  analysis,  and  those  that  may  result  in 
intractable  run  times.  Again,  it  is  reiterated  that  the  end  goal  is  performing 
production  computing  involving  more  realistic  applications. 


Table  9.  FPC  computational  cost  comparison,  GSA  configuration. 


FPC 

Mesh 

Resolution 

(cm) 

Approximate 
No.  of  Cells 
(millions) 

Anaiysis 

End  Time 
(msec) 

No.  of 
Processors 

Waii  Ciock 
Time  (hr) 

Totai 

Processor- 

Hours 

CTH 

0.5 

395 

60.5 

1024 

147 

150,528 

DYSMAS 

1.0 

136 

62.7 

1016 

161 

16,256 

Loci/BLAST 

2.5 

52.4 

70 

200 

156.5 

31,300 

RAGE 

0.5  (AMR) 

13  -  279 

52.6 

128  - 1024 

429 

282,2802 

SHAMRC 

(Coarse) 

1.0 

256 

110 

64 

31 

1,984 

SHAMRC 

(Medium-HS) 

0.5 

1024 

110 

128  -  256 

93.5 

20,992 

1  Reported  as  CPU  time. 

2  Reported  in  Processor-CPU  hours. 


Table  10.  FPC  computational  cost  comparison, 
8-ftx8-ft  C2SQ  configuration. 


FPC 

Mesh 

Resoiution 

(cm) 

Approximate 
No.  of  Ceiis 
(miiiions) 

Anaiysis 
End  Time 
(msec) 

No.  of 
Processors 

Waii 

Ciock 

Time  (hr) 

Totai 

Processor- 

Hours 

CTH 

0.5 

473 

67.4 

1024 

193 

197,632 

DYSMAS 

1.0 

166 

150 

1016 

31.91 

32,410 

Loci/BLAST 

2.5 

73.1 

70 

200 

194.9 

38,980 

RAGE 

0.5  (AMR) 

39  -  509 

57.2 

256  - 1536 

653 

611,4902 

SHAMRC 

(Coarse) 

1.0 

272 

110 

256 

26.6 

6,810 

SHAMRC 

(Medium-HS) 

0.5 

1088 

110 

256 

68.9 

17,638 

1  Reported  as  CPU  time. 

2  Reported  in  Processor-CPU  hours. 
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Table  11.  FPC  computational  cost  comparison, 
8-ftx8-ft  C2SQ  configuration  with  box. 


FPC 

Mesh 

Resolution 

(cm) 

Approximate 
No.  of  Cells 
(millions) 

Anaiysis 

End 

Time 

(msec) 

No.  of 
Processors 

Waii  Ciock 
Time  (hr) 

Totai 

Processor- 

Hours 

CTH 

0.5 

473 

70 

1024 

199 

203,776 

DYSMAS 

1.0 

166 

150 

1016 

37.51 

38,100 

Loci/BLAST 

2.5 

69.9 

70 

200 

240 

48,000 

RAGE 

0.5  (AMR) 

39  -  478 

60 

256  -  2048 

643 

661,9802 

SHAMRC 

(Coarse) 

1.0 

272 

110 

128 

52 

6,656 

SHAMRC 

(Medium-HS) 

0.5 

1088 

110 

128  -  256 

101.3 

19,808 

1  Reported  as  CPU  time. 

2  Reported  in  Processor-CPU  hours. 


There  is  one  last  point  to  make  regarding  the  SHAMRC  calculations.  Mul¬ 
tiple  sets  of  calculations  were  conducted  for  each  BLS  configuration.  Two 
sets,  denoted  coarse  and  medium  with  half-symmetry  (HS),  were 
conducted  for  the  GSA  and  8-ftx8-ft  C2SQ  configurations.  Two  additional 
sets  of  calculations  were  conducted  for  the  third  case  (medium  fully  3-D 
and  fine  half-symmetry).  Only  those  calculations  that  can  be  compared 
across  all  three  BLS  configurations  are  included  in  the  tables.  The  reader  is 
directed  to  Appendix  E  for  details  on  the  other  SHAMRC  calculations. 

One  would  like  to  use  a  metric  like  the  grind  time  to  compare  expense  on  a 
cycle-by-cycle  basis  (Attaway  et  al.  1998).  The  grind  time  tgrind  is  defined  as 

^  _  ^execution  //i  \ 

^arind-  N,, us  N cycles 

where  texecution  is  the  execution  time,  Nceiis  is  the  total  number  of  cells  in  the 
problem,  and  Ncycies  is  the  total  number  of  cycles  (or  time-steps)  taken  in 
the  calculation.  It  provides  a  measure  of  the  CPU  time  required  to  perform 
a  calculation  for  a  single  cell  over  one  time-step.  The  grind  time  is  problem 
dependent,  and  in  the  case  of  transient  problems,  it  can  be  affected  by 
variable  time-stepping.  When  comparing  grind  times  for  different  codes,  it 
is  necessary  that  the  problem  setups  be  identical  (or  relatively  similar).  For 
this  work,  that  means  a  comparable  physical  space  encompassing  the  BLS 
in  each  calculation  as  well  as  identical  approaches  for  modeling  interior 
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features  such  as  the  striker  and  grill.  This  is  not  the  case  for  the  BLS  mod¬ 
eling  effort.  Further,  code-to-code  comparisons  based  on  grind  time  (or 
any  other  computational  efficiency  measure)  should  all  be  performed  on 
the  same  HPC  platform.  This  is  not  the  case  for  the  BLS  modeling  effort 
either.  Thus,  it  is  necessary  to  develop  an  ad  hoc  approach  to  obtain  cost 
comparisons  on  as  consistent  a  basis  as  possible. 

After  reviewing  the  data  in  the  tables,  the  decision  was  made  to  draw 
comparisons  on  a  processor-hour  (PH)  basis  scaled  to  an  end  analysis  time 
of  100  msec.  This  measure  is  defined  as 

PHioo  msec  =  P  (100  msec)  (2) 

where  the  final  value  of  PH  and  end  analysis  time  (tend)  are  specified  in 
Tables  8  through  lo.  The  scale  factor  beta  (|3)  is  applied  based  on  symme¬ 
try  conditions.  For  a  fully  3-D  calculation,  beta  is  set  to  1.  For  a  half¬ 
symmetry  calculation,  beta  is  set  to  2.  The  goal  is  to  equate  all  PH 
estimates  for  a  fully  3-D  calculation.  It  is  assumed  that  the  problem  size 
will  grow,  and  that  the  number  of  processors  required  capturing  the 
increased  problem  size  doubles.  The  scaled  costs  are  outlined  in  Table  12. 
Note  the  PH  value  is  expressed  in  terms  of  processor-kilohours  (per 
1,000  hr).  This  was  done  to  simplify  comparisons  and  is  in  line  with  the 
goal  of  an  order  of  magnitude  estimate.  The  reader  should  not  expect 
anything  more  than  a  rough  estimate,  given  the  disparate  setups  taken  by 
each  of  the  modelers. 


Table  12.  FPC  costs  scaled  to  100  msec  and  fully  3-D  calculation. 


FPC 

Processor  -  K  hours  to  reach  100  msec 

GSA 

Configuration 

8-ftx8-ftC2SQ 

Configuration 

8-ftx8-ftC2SQ 
Configuration 
with  Box 

CTH 

498 

586 

582 

DYSMAS 

26 

22 

25 

Loci/BLAST 

45 

56 

69 

RAGE 

1075 

2138 

2207 

SHAMRC  (Coarse) 

2 

7 

6 

SHAMRC  (Medium-HS) 

42 

33 

36 
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It  is  apparent  that  SHAMRC,  DYSMAS,  and  Loci/BLAST  are  all  cost  com¬ 
petitive  and  fall  into  the  first  tier.  Both  DYSMAS  and  the  SHAMRC 
(coarse)  calculations  utilized  a  i.o-cm  resolution  mesh.  On  this  basis,  one 
might  rank  SHAMRC  as  the  best  cost  performer  for  the  codes  evaluated. 
The  cost  for  CTH  is  an  order  of  magnitude  greater  than  any  of  the  first  tier 
codes.  Both  CTH  and  the  SHAMRC  (medium)  calculations  utilized  a 
0.5-cm  resolution  mesh,  with  CTH  requiring  about  12  to  16  times  more 
resources.  Recall  that  the  PH  estimate  is  a  combination  of  time  and 
processors,  with  the  combination  of  the  two  representing  the 
computational  resources  needed  to  conduct  the  calculation. 

The  RAGE  calculations  were  exceedingly  expensive.  This  is  in  line  with  our 
experience  modeling  explosive  detonations  using  the  SAIC  Adaptive  Grid 
Generation  (SAGE)  code  at  the  ERDC  (Archer  et  al.  2005).  The  AMR 
capability  in  SAGE  performs  exceedingly  well  in  capturing  explosive  deto¬ 
nation  processes  as  well  as  the  near  field  blast  behavior;  however,  run 
times  increase  substantially  when  the  physical  space  being  modeled  is 
large.  Examples  involving  a  large  physical  space  include  both  open-air  and 
internal  detonations.  SAGE  (and  RAGE)  has  the  capability  to  map  a  i-DS 
or  2-DC  calculation  into  a  3-D  problem  setup,  thereby  reducing  the  mesh 
resolution  requirements  for  the  latter.  Unfortunately,  this  capability  could 
not  be  utilized  in  the  BLS  modeling  that  required  that  all  RAGE 
calculations  be  3-D  from  the  start. 
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4  Summary  and  Recommendations 

An  integrated  experimental  and  computational  program  was  conducted  to 
evaluate  several  FPCs  for  modeling  airblast  environments  typical  of  that 
encountered  in  the  ERDC  BLS.  The  FPCs  considered  were  CTH,  DYSMAS, 
Loci/BLAST,  RAGE,  and  SHAMRC.  These  codes  afforded  a  range  of 
numerical  solvers  and  EOS  options  for  the  modeling  effort.  The  FPCs  were 
evaluated  against  data  generated  in  the  BLS.  Three  experimental 
configurations  were  considered  for  the  validation.  These  included  two 
empty  configurations,  with  comparisons  drawn  against  reflected  pressures 
measured  at  the  downstream  target  wall.  The  third  included  a  hox-like 
structure  in  the  flow  field.  The  structure  was  instrumented  on  the  exposed 
sides  to  provide  a  more  challenging  scenario  for  the  model  comparisons. 
The  purpose  of  the  integrated  program  was  to  assess  computational  accu¬ 
racy  and  cost  of  the  aforementioned  FPCs,  as  well  as  identify  any  short¬ 
comings  in  the  physics  modeling  and  areas  for  future  improvement. 

A  qualitative  assessment  of  accuracy  was  made  based  on  visual  compari¬ 
sons  between  the  calculated  and  measured  waveforms  for  the  three  BLS 
configurations  modeled.  Good  overall  correlation  with  the  measured  data 
was  noted  for  CTH,  DYSMAS,  Loci/BLAST,  and  SHAMRC.  The  relative 
ranking  of  these  codes  based  on  accuracy  is  problem  dependent,  but  cer¬ 
tainly  results  from  any  of  these  would  be  considered  reasonable.  A  sepa¬ 
rate  assessment  was  made  based  on  computational  cost.  Run  times  for 
DYSMAS,  Loci/BLAST,  and  SHAMRC  were  quite  acceptable  for  produc¬ 
tion  computing.  CTH  was  found  to  be  very  expensive,  with  the  computa¬ 
tional  cost  being  an  order  of  magnitude  greater  than  that  of  the  other  three 
codes. 

RAGE  was  found  to  perform  poorly  overall  (in  both  accuracy  and  cost).  It 
is  conjectured  that  there  is  an  issue  with  the  ideal  gas  EOS  implementation 
that  should  be  addressed.  This  is  discussed  in  Appendix  D.  The  Leidos 
version  of  RAGE  was  utilized  in  this  work,  and  it  is  recommended  that  a 
follow-on  analysis  be  conducted  with  the  LANL  version.  This  would  be 
helpful  to  diagnose  the  issue  and  possibly  pinpoint  the  source  of  the 
problem.  The  follow-on  work  could  be  very  limited  in  scope.  A  simple  2- 
DC  calculation  would  suffice  and  would  provide  a  quick  assessment  of 
whether  the  issue  persists  in  the  LANL  version. 
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All  of  the  calculations  conducted  for  this  report  assumed  inviscid  flow. 
Consequently,  physical  mechanisms  such  as  turbulence  and  boundary 
layer  effects  were  not  included  in  the  simulations.  Review  of  the  user 
documentation  suggests  that  both  Loci/CHEM,  a  forerunner  of  Loci/ 
BLAST  that  has  an  implicit  solver  built  under  the  same  framework,  and 
SHAMRC  have  the  capability  to  handle  viscosity  and/or  turbulence. 
Although  the  inviscid  flow  calculations  considered  here  suggest  that  vis¬ 
cosity  has  a  second  order  effect,  this  may  not  always  be  the  case  as  the 
experimental  program  moves  forward.  It  is  recommended  that  a  follow-on 
analysis  that  includes  viscous  effects  be  conducted  with  one  of  the  BLS 
configurations  considered  here  to  examine  if  the  solutions  change.  This 
would  provide  an  initial  validation  test  of  the  viscous  modeling  and  pro¬ 
vide  insight  into  the  computational  costs. 

The  BLS  testing  considered  in  the  study  was  based  on  a  single  set  of  driver 
conditions.  It  is  recommended  that  extended  testing  be  conducted 
involving  the  third  BLS  configuration  (with  box).  This  testing  could 
maintain  the  same  configuration  with  variations  in  the  driver  conditions. 
This  would  allow  the  accuracy  of  the  models  to  be  evaluated  over  a  range 
of  loading  conditions.  It  is  further  recommended  that  the  testing  be 
extended  to  consider  a  two-box  configuration.  Here  the  emphasis  would  be 
on  capturing  the  diffraction  and  vortex  shedding  noted  about  the  first  box 
and  subsequent  pressure  loading  on  the  face  of  the  second  box.  The  spac¬ 
ing  between  boxes  would  influence  the  degree  of  shock  interaction  and 
subsequently  the  loads  induced  on  the  box  faces.  A  methodical  approach 
would  be  to  consider  two  boxes  of  the  same  height  and  vary  the  spacing  in- 
between  for  a  fixed  set  of  driver  conditions.  The  driver  conditions  could 
then  be  varied  to  induce  a  range  of  incident  shocks  conditions.  This  would 
allow  for  the  codes  to  be  evaluated  over  a  range  of  conditions  using  a  rela¬ 
tively  simple  problem  setup. 
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Appendix  A:  CTH  Modeling 

CTH  Background 

Comprehensive  descriptions  of  the  CTH  code  are  provided  in  Crawford 
et  al.  (2013b),  Bessette  and  Britt  (2011),  and  McGlaun  et  al.  (1990).  The 
following  is  a  summary  of  those  descriptions.  CTH  is  a  multi-material, 
Eulerian,  large  deformation,  strong  shock  wave,  solid  mechanics  code 
developed  at  Sandia  National  Laboratories.  It  has  been  in  wide  use  within 
the  DoD  community  for  several  decades.  CTH  falls  into  a  category  of  struc¬ 
tured  mesh  codes  in  which  the  reference  mesh  is  fixed  in  space  and  aligned 
with  the  major  coordinate  axes.  Three-dimensional  (3-D)  rectangular 
meshes,  two-dimensional  rectangular  (2-DR)  and  cylindrical  (2-DC) 
meshes,  and  one-dimensional  rectilinear  (i-DR),  cylindrical  (i-DC),  and 
spherical  (i-DS)  meshes  are  available  to  the  user.  CTH  also  has  adaptive 
mesh  refinement.  A  3-D  mesh  is  composed  of  hexahedral  cells  that  can 
vary  in  size  and  aspect  ratio  along  the  major  coordinate  axes.  As  described 
by  Bessette  and  Britt  (2011),  CTH  utilizes  a  two-step  solution  procedure 
for  the  conservation  equations.  It  first  takes  a  Lagrangian  step,  which 
allows  the  reference  mesh  to  deform.  The  equations  of  motion  are 
integrated  forward  in  time  with  the  material  response  updated  subject  to  a 
constitutive  model.  The  integration  is  performed  explicitly  with  individual 
time-steps  subject  to  a  Courant-Friedrichs-Lewy  (CFL)  condition  on 
stability,  as  well  as  limits  on  volume  transport  between  cells.  The  CTH 
default  CFL  values  of  0.55  for  3-D  and  0.6  for  i-D  and  2-D  were  used  for 
all  calculations. 

The  Lagrangian  step  is  followed  by  a  remap  or  advection  step  that  maps 
data  (i.e.,  the  volume  flux,  mass,  momentum,  energy,  and  material  state 
data)  from  the  deformed  Lagrangian  configuration  back  into  the  fixed 
Eulerian  frame  of  reference.  In  this  step,  a  donor-acceptor  relationship  is 
defined  where  a  donor  cell  is  said  to  advect  material  to  an  acceptor  cell. 
First,  the  volume  flux  from  the  donor  to  an  acceptor  cell  is  determined. 
Then  the  material  mass  and  internal  energy  are  advected.  Finally,  the 
momentum  and  kinetic  energy  are  advected  to  the  reference  mesh.  Advec¬ 
tion  is  determined  using  a  second-order  accurate  scheme  attributed  to  Van 
Leer  (1997).  An  interface-tracking  algorithm  is  used  to  reconstruct  the 
locations  of  material  interfaces  (as  well  as  free  surfaces)  based  on  the 
volume  fraction  of  materials  in  the  upstream,  downstream  (acceptor),  and 
donor  cells. 
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CTH  supports  a  wide  range  of  constitutive  models  appropriate  for  strong 
shock,  large  deformation  calculations.  It  has  models  for  multiphase,  elas¬ 
tic,  viscoplastic,  porous,  and  explosive  materials,  with  model  coefficients 
provided  in  a  library  supplied  with  the  code.  The  stress  response  in  a 
material  is  decomposed  into  its  dilatational  and  deviatoric  behaviors.  The 
dilatational  (or  volumetric)  response  is  described  using  an  EOS,  while  a 
strength  model  is  used  to  describe  the  deviatoric  (or  distortional) 
response.  In  the  high-pressure  regime,  SESAME  tabular  and  analytic  EOS 
models  are  used  to  replicate  the  nonlinear  behavior  of  materials.  The 
SESAME  table  (Hertel  and  Kerley  2006)  is  intended  to  cover  a  broad 
density-temperature  range  so  that  it  can  be  used  in  many  different 
applications,  e.g.,  to  treat  high  compressions,  large  expansions,  and  high 
temperatures.  SESAME  can  model  solid,  liquid,  vapor,  liquid-vapor,  solid- 
liquid,  and  solid-solid  phase  changes.  Hertel  and  Kerley  (2006)  stated  that 
the  SESAME  table  “allows  the  use  of  sophisticated  models  that  are  too 
complicated  to  be  incorporated  into  analytic  formulas.  A  good  tabular  EOS 
gives  valid  results  over  a  much  wider  density-temperature  range  than  the 
analytic  EOS  models.” 

High-explosive  detonation  can  be  modeled  using  programmed  burn  or 
several  reactive  burn  models,  e.g.,  the  History  Variable  Reactive  Burn 
(HVRB)  model  (Hertel  and  Kerley  2006)  developed  at  Sandia.  With  the 
program  burn  option,  “an  ideal  detonation  wave  front  is  made  to  propa¬ 
gate  at  a  given  velocity  from  initiation  points  specified  by  the  user”  (Hertel 
and  Kerley  2006).  In  contrast  to  programmed-burn  models,  the  material 
decomposition  is  included  in  the  EOS  of  reactive  burn  models,  and  time 
evolution  of  the  reaction  is  described  by  a  rate  equation  (Hertel  and  Kerley 
2006).  Either  the  Jones-Wilkins-Lee  (JWL)  analytic  (Hertel  and  Kerley 
2006)  or  the  SESAME  tabular  equations  of  state  can  model  the  high- 
explosive  reaction  products  (Crawford  et  al.  2013b). 

In  the  problem  development,  the  user  may  define  material  or  part  geome¬ 
tries  using  simple  primitive  shapes,  which  can  be  combined  (or  sub¬ 
tracted)  as  needed  to  develop  a  final  geometric  entity.  As  an  alternative, 
computer  aided  design  (CAD)  geometry  may  be  inserted  using  a  limited 
number  of  file  formats  (Crawford  et  al.  2013b),  i.e.,  ASCII  versions  of  a 
Pro/Engineer  tetrahedral  or  shell  mesh,  an  Exodus  tetrahedral  or 
hexahedral  mesh,  or  an  STL  file.  Boundary  conditions  (BCs)  are  applied  to 
the  borders  of  the  fixed  reference  mesh.  CTH  supports  several  BC  types, 
including  reflective  or  symmetry,  transmissive,  and  outflow  (Crawford 
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et  al.  2013b).  Transmissive  BCs  are  sound  speed  absorbing  and  are  typi¬ 
cally  used  to  model  semi-infinite  domains.  However,  they  can  be  problem¬ 
atic  for  strong  shock  interactions  where  wave  speeds  exceed  the  acoustic 
velocity. 

In  CTH,  artificial  bulk  viscosity  terms  (linear  and  quadratic)  assist  with  the 
handling  of  strong  discontinuities  in  stress,  such  as  those  found  in  a  shock 
front.  The  quadratic  term  is  used  to  smear  a  shock  front  across  several 
cells.  The  linear  term  is  used  to  minimize  excess  noise  (“ringing”)  in  the 
simulation.  A  linear  term  that  is  too  large  can  cause  excessive  smearing  of 
the  discontinuity  (Crawford  et  al.  2013b). 

Within  the  CTH  family  of  codes,  a  rezoning  program  (Yarrington  2006) 
allows  an  analyst  to  modify  a  simulation  at  a  specified  time.  The  user  may 
change  the  zoning,  orientation,  and  the  number  of  dimensions.  For  exam¬ 
ple,  a  i-DS  calculation  of  an  aboveground  detonation  with  0.2-cm  cell  sizes 
could  be  run  out  until  the  shock  front  nears  the  ground  surface.  Those 
results  can  be  rezoned  into  a  3-D  mesh  with  nominal  i-in.  cells,  and  the 
calculation  continued  out  to  later  times.  This  process  allows  the  analyst  to 
use  a  much  finer  mesh  for  the  detonation  region  and  a  more  affordable 
mesh  size  for  the  subsequent  3-D  calculation,  resulting  in  a  great  savings 
of  time  and  computer  resources  versus  running  the  entire  simulation  in  3- 
D.  Rezoning  was  not  used  in  this  project,  due  to  past  problems  with  the 
software. 

CTH  requires  all  input  in  CGS  units,  i.e.,  centimeter  (cm),  gram  (g),  sec¬ 
ond  (sec),  and  electron  volts  for  temperature  (one  eV  =  11,600  deg  Kelvin). 
Pressures  in  this  unit  system  are  in  dynes/cm^. 

Analysis  overview 

All  of  the  CTH  calculations  were  run  on  the  Cray  XE6  located  at  the  ERDC 
DoD  Supercomputing  Resource  Center,  commonly  referred  to  as  Garnet. 
The  Cray  XE6  has  4,716  compute  nodes  with  32  cores  per  node.  Each 
compute  node  has  a  core  speed  of  2.5  GHz  and  64  GB  of  accessible 
memory. 

A  stepwise  approach  was  taken  to  evaluate  the  shock  environments  for  this 
project.  This  involved  a  series  of  i-DS,  2-DC,  and  eventually  3- 
D  calculations,  with  a  progressively  increasing  level  of  difficulty.  Both  the 
SESAME  and  ideal  gas  EOSs  were  used  to  approximate  the  properties  of 
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the  air.  The  walls  of  the  BLS  were  modeled  as  either  rigid  or  deformable 
materials.  Deformable  walls  were  assigned  the  properties  of  lead  due  to  its 
low  wave  speed  and  the  desire  to  achieve  a  time-step  as  large  as  possible. 
The  lead  properties  were  approximated  with  the  Mie-Griineisen  EOS  and 
with  an  enhanced  von  Mises  yield  strength.  The  elastic  parameters  were 
not  altered.  The  striker  was  modeled  in  all  of  the  3-D  calculations  and  in 
six  of  the  fifteen  2-DC  simulations.  Other  than  its  nominal  3-in.  width,  the 
grill  was  neither  modeled  in  the  2-DC  simulations,  due  to  the 
axisymmetric  geometry,  nor  in  the  3-D  calculations,  due  to  the  fine  mesh 
resolution  required  to  accurately  capture  the  individual  horizontal  and 
vertical  i/8-in.-thick  members. 

CTH  requires  input  for  time-step  control,  restarts,  cell  thermodynamics, 
convection,  discards,  tracer  points,  and  spatial  output.  Discards  are  used 
to  improve  the  time-step  and/or  prevent  a  calculation  from  stopping  due 
to  “material  that  is  acquiring  non-physical  states.”  Spatial  output  (e.g., 
pressure,  density,  or  temperature  plotted  over  the  computational  domain) 
was  saved  to  Spymaster  files  (Crawford  et  al.  2013a).  These  files  were  post- 
processed  with  either  CTH’s  SPYPLT  software  or  with  the  open-source 
visualization  application  ParaView  (Kitware  2015). 

Time-history  output  locations  in  CTH  are  defined  with  tracer  input.  Each 
tracer  may  be  fixed  in  space  or  allowed  to  move.  For  the  BLS  simulations, 
all  of  the  tracers  were  fixed.  As  explained  in  Thompson  and  Kmetyk 
(1999),  plot  variables  consist  of  either  point  (e.g.,  pressure,  density, 
temperature),  material  (e.g.,  energy,  momentum),  or  global  (e.g.,  time- 
step,  cycle  number,  mass)  information.  Approximately  70  variables  are 
available  for  post-processing.  Care  must  be  taken  to  avoid  placing  tracers 
at  locations  where  mixed  cells  (more  than  one  material  per  cell)  may 
develop.  A  mixed  cell  with  air  and  solid  will  typically  produce  undesired 
output  for  the  user. 

Early  1-DS  calculations 

The  i-DS  calculations  were  conducted  for  two  purposes:  (1)  to  evaluate 
mesh  refinement,  and  (2)  to  evaluate  the  two  EOS  options  for  air.  Both  the 
SESAME  and  ideal  gas  EOSs  were  used  to  represent  the  air.  The  i-DS 
calculations  were  very  cost-effective.  For  example,  a  typical  run  with 
0.2-cm  resolution,  a  maximum  range  of  20  ft,  with  approximately 
3,100  cells,  and  a  simulation  time  of  20  m/s  required  less  than  40  sec  of 
wall-clock  time  using  32  cores  (or  one  compute  node). 
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Figure  A  shows  the  difference  in  calculated  peak  overpressure  as  a  func¬ 
tion  of  range  for  meshes  having  cell  sizes  of  0.2  and  2.0  cm.  The  source 
region  in  the  2.0-cm  resolution  calculation  had  a  finer  mesh  (0.2-cm  cell 
size)  to  ensure  consistency  in  the  comparisons.  At  the  6.6-ft  range,  the 
peak  overpressure  in  the  finer  mesh  was  6  percent  greater  than  the  coarser 
one.  At  the  3.8-ft  range,  the  increase  was  a  7.5  percent.  In  Figure  A2,  the 
peak  overpressure  versus  range  is  shown  for  calculations  using  the 
SESAME  and  ideal  gas  EOSs  for  air.  At  the  6.6-ft  range,  the  peak  pressure 
associated  with  the  ideal  gas  EOS  was  10.9  percent  greater  and  at  the 
3.8-ft  range  it  was  10.7  percent  greater. 


Figure  Al.  Results  of  1-DS  CTH  calculations 
comparing  mesh  resolution. 
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Figure  A2.  Results  of  1-DS  CTH  calculations 
comparing  EOS  options. 
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Early  2-DC  calculations 

The  2-DC  CTH  calculations  of  the  GSA  configuration  (Case  i)  were  con¬ 
ducted  to  evaluate  (i)  the  effect  of  mesh  refinement  on  peak  overpres¬ 
sures,  (2)  different  air  EOSs,  (3)  the  effect  of  the  striker  on  downstream 
pressures,  and  (4)  the  planarity  of  the  shock  in  an  axisymmetric  geometry. 
The  2-DC  calculations  used  nominal  0.1-,  0.2-,  or  0.5-cm  cell  sizes.  Reflec¬ 
tive  BCs  were  applied  to  the  axis  of  symmetry  and  transmissive  BCs  were 
applied  to  the  three  borders  of  the  computational  domain.  The  walls  of  the 
BLS  were  modeled  as  rigid  materials.  The  air  was  modeled  with  either  the 
SESAME  or  ideal  gas  EOS.  The  driver  air  was  assigned  an  initial  pressure 
of  1,298  psig  and  an  initial  temperature  of  77°F,  and  the  ambient  air  initial 
values  of  14.7  psia  and  77OF. 

Figures  A3  through  A8  show  the  results  from  a  nominal  0.2-cm  resolution 
mesh  and  the  SESAME  EOS  for  air.  Pressure  state  plots  at  12  msec 
(a)  without  and  (b)  with  the  striker  are  shown  in  Figure  A3.  The  pressure 
environment  between  the  driver  and  upstream  edge  of  the  CVC-2  section 
has  been  significantly  perturbed  by  the  presence  of  the  striker.  However, 
the  environments  from  the  upstream  edge  of  the  cone  section  to  the  right 
edge  of  the  image  are  qualitatively  similar.  The  expansion  of  the  air  in  the 
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cone  section  appears  to  mitigate  the  heterogeneous  state  of  the  pressure 
field  induced  by  the  presence  of  the  striker. 

Figure  A4  shows  the  early-time  pressure  environment  at  3,  5,  and  7  msec 
with  the  striker  installed.  At  3  msec,  the  shock  has  reached  the  entrance  of 
the  cone  section  with  a  non-planar  front,  and  high  pressure  air  is  venting 
out  of  the  system.  At  5  msec,  the  initial  shock  is  approximately  one-quarter 
down  the  length  of  the  cone  section  and  is  still  non-planar.  Within  the  CVC 
section,  the  pressure  environment  is  very  heterogeneous  and  shows 
multiple  reflected  shocks.  The  initial  shock  at  7  msec  is  approximately  at 
the  center  of  the  cone  section  and  still  exhibits  some  curvature.  The  com¬ 
plex  nature  of  the  environment  is  still  present  within  the  CVC  section. 


Figure  A3.  CTH  2-DC  pressure  state  plots  at 
12  msec,  with  and  without  striker. 


(a)  Without  striker 


{b)  With  striker 
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Figure  A4.  CTH  2-DC  pressure  state  plots  at  varied  times. 


(a)  3  msec 


(b)  5  msec 


(c)  7  msec 
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Pressure  and  velocity  magnitude  are  shown  at  12  msec  in  Figure  A5(a)  and 
(b),  respectively.  Each  image  shows  two  views;  on  the  left  the  entire  BLS 
and  on  the  right  the  upstream  end  of  the  BLS.  The  initial  shock  front  has 
entered  the  TRi  and  TR2  sections  of  the  BLS  at  this  time.  The  pressure 
environment  behind  the  front  is  more  homogeneous  than  at  the  earlier 
times.  However,  the  velocity  magnitude  plot  exhibits  a  more  complex 
structure  than  the  pressure.  High  pressure  air  escaping  from  the  vents  is 
clearly  evident  in  this  image. 


Figure  A5.  CTH  2-DC  state  plots  at  12  msec. 


(a)  Pressure  at  12  msec 


(b)  Velocity  magnitude  at  12  msec 


The  pressure  environment  in  the  BLS  at  15  and  27  msec  is  captured  in 
Figure  A6.  At  15  msec,  the  initial  shock  front  is  planar,  and  a  mass  of  air  at 
nominal  constant  pressure  follows  the  front.  At  20  msec  (no  image  pro¬ 
vided),  the  initial  shock  front  reached  the  bulkhead  at  the  end  of  the  BLS. 
At  27  msec,  the  reflected  shock  front  is  positioned  approximately  at  the 
middle  of  the  TR-Cascade  section.  The  reflected  front  at  this  time  is  non- 


ERDC/GSL  TR-16-31 


57 


planar  due  to  the  interaction  of  the  front  with  the  mass  of  air  still  moving 
downstream.  This  image  also  shows  intermittent  focusing  of  higher  pres¬ 
sure  air  along  the  centerline  of  the  BLS. 


Figure  A6.  CTH  2-DC  pressure  state  plots  at  later  times. 


(a)  Pressure  at  15  msec 


Figure  A7  shows  time-history  traces  at  gage  locations  P3,  P4,  and  PE. 
Image  (a)  shows  a  comparison  plot  of  the  traces  from  the  calculation.  The 
P3  trace  contains  significant  noise  compared  to  the  other  two.  This  noise  is 
attributed  to  the  tracer  location  being  located  in  a  cell  2.25  in.  away  from 
the  axis  of  symmetry.  The  CTH  time-history  traces  at  gage  locations  P4 
and  PE  are  compared  to  the  experimental  records  in  images  (b)  and  (c), 
respectively.  The  CTH  data  were  shifted  -4.59  msec  to  match  the  time  of 
arrival  of  the  experimental  records.  The  peak  pressures  in  the  two  CTH 
traces  were  greater  than  the  measured  and  the  overall  impulse  was  not 
captured  well  in  CTH. 
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Figure  A7.  Gauge  comparisons  for  2-DC  CTH  caicuiation 
of  GSA  configuration. 


(a)  Caicuiated  oniy 


(b)  Gauge  P4 


(c)  Gauge  PE 


Overview  of  the  3-D  calculations 

In  general,  the  3-D  CTH  calculations  of  the  BLS  were  conducted  using 
half-symmetry,  a  nominal  0.5-cm  cell  size,  the  SESAME  EOS  for  air,  and 
deformable  walls.  There  was  one  exception  that  utilized  the  same 
configuration  and  resolution;  however,  the  air  was  modeled  as  an  ideal 
gas.  As  described  in  the  main  text,  three  configurations  were  simulated, 
the  GSA  configuration,  the  8-ftx8-ft  C2SQ  configuration,  and  the  8-ftx8-ft 
C2SQ  configuration  with  box.  Each  subsequent  configuration  added  more 
complexity  to  the  modeling  effort. 

Reflective  BCs  were  applied  to  the  axis  of  symmetry,  and  transmissive  BCs 
were  applied  to  the  three  borders  of  the  computational  domain.  The  BLS 


Impulse  (psi-msec) 


ERDC/GSL  TR-16-31 


59 


geometry  was  created  in  Cubit  (a  general  pre-processing  tool;  Sandia 
National  Laboratories  2015),  and  the  major  components  output  to  individ¬ 
ual  Exodus  files.  CTH  imports  the  Exodus  files  and  extracts  the  volume  of 
each  component.  It  then  simply  superimposes  the  prescribed  computa¬ 
tional  mesh  on  top  of  the  volumes.  Eigure  A8  shows  images  of  the  compo¬ 
nents  for  the  8-ftx8-ft  C2SQ  configuration  with  box.  The  four  colors 
represent  the  four  volumes  imported  into  CTH.  A  close-up  of  the  as- 
modeled  striker  geometry  is  shown  in  Eigure  A9. 

GSA  configurations,  Tests  30a  and  30b 

The  simulation  for  the  GSA  configuration  used  395  million  cells  and 
required  147  hr  of  wall-clock  time  using  1,024  cores.  No  restarts  were 
required.  The  resulting  pressure  environments  from  the  3-D  calculation 
were  not  significantly  different  from  the  previously  described  2-DC 
results;  thus,  only  the  CTH  time-history  traces  for  gages  P3,  P4,  and  PE 
are  presented  herein.  These  traces  were  shifted  in  time  -4.59  msec  and  are 
compared  to  the  measured  data  in  Eigure  Aio.  In  general,  there  was  good 
agreement  between  the  CTH  pressure  and  impulse  traces  and  the 
experimental  results.  The  one  exception  was  gage  P3  (see  Figure  Aio(a)) 
where  the  CTH  pressures  after  40  msec  exhibited  two  excursions  above 
the  experimental  records,  which  in  turn  produced  a  higher  impulse  trace 
after  40  msec.  The  CTH  pressure  traces  also  indicate  that  the  numerical 
shock  front  was  planar  when  it  reached  the  bulkhead,  as  all  of  the  traces 
had  the  same  time  of  arrival.  Figure  All  compares  the  2-DC  (0.2-cm 
resolution)  and  the  3-D  (0.5-cm  resolution)  pressures  and  impulses  for  the 
P4  and  PE  gage  locations.  The  oscillations  in  the  2-DC  traces  are 
attributed  to  the  finer  mesh. 
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Figure  A8.  Geometry  for  the  8-ftx8-ft  C2SQ 
configuration  with  box. 


(a)  image  of  half-symmetry  model  used  in  CTH 


(b)  Quarter-symmetry  image  of  model 


(d)  Quarter-symmetry  image  of  instrumented  target  box 
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Figure  A9.  Three  views  of  as-modeied  striker  geometry. 


(a)  Downstream  view  (b)  Upstream  view 


(c)  Side  view 


Overpressure  (psi)  Overpressure  (psi) 
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Figure  AlO.  Gauge  comparisons  for  CTH,  GSA  configuration. 


(a)  Gauge  P3  (b)  Gauge  P4 


Time  (msec) 


(c)  Gauge  PE 


Figure  All.  Comparisons  from  2-DC  and  3-D  CTH  caicuiations. 


Time  (msec)  Time  (msec) 

(a)  Gauge  P4  (b)  Gauge  PE 


Impulse  (psi-msec)  Impulse  (psi-msec) 
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Calibration  tests  with  8-ftx8-ft  C2SQ  configuration 

The  CTH  calculation  of  the  8-ftx8-ft  C2SQ  configuration  used  473  million 
cells  and  required  193  hr  of  wall-clock  time  using  1,024  cores.  One  restart 
was  required  to  reach  72  msec  of  simulation  time.  Since  this  simulation 
used  half  symmetry,  gauge  locations  on  the  opposing  half  were  mirrored 
about  the  plane  of  symmetry.  The  CTH  pressure  and  impulse  traces  for  the 
target  wall  gauges  are  compared  with  the  measured  data  in  Figure  A12. 

The  CTH  time-histories  were  shifted  -4.80  msec  to  match  the  measured 
arrival  time. 

With  the  exception  of  gauges  P7  and  P8,  there  was  excellent  agreement 
between  the  calculations  and  measured  data.  The  discrepancy  for  gauge  P7 
is  large;  however,  the  measured  pressure-time-history  at  this  particular 
gauge  differs  markedly  from  that  measured  at  the  other  gauge  locations. 
Although  the  measurement  appears  repeatable  from  test-to-test,  one 
might  still  consider  the  data  for  gauge  P7  suspect  because  it  is  so  out-of- 
line  with  the  other  measurements. 

The  CTH  pressure  traces  indicate  that  the  numerical  shock  front  was 
planar  when  it  reached  the  target  wall,  as  all  of  the  traces  had  the  same 
nominal  TOA.  Figure  A13  presents  the  CTH  and  experimental  traces  for 
the  side-on  overpressure  gauge  located  in  the  TR2  section  of  the  BLS.  The 
initial  pressure  traces  agree  between  the  TOA  and  17  msec,  after  which  the 
CTH  trace  falls  below  the  experimental  records.  At  26  msec,  the  CTH  trace 
follows  the  experimental  data  again.  CTH  overpredicts  the  second  peak  at 
40  msec. 


Overpressure  (psi)  Overpressure  (psi)  Overpressure  (psi) 
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Figure  A12.  Gauge  comparisons  for  CTH,  8-ftx8-ft 
C2SQ  configuration. 


Time  (msec) 

(a)  Gauge  PI 


Time  (msec) 

(c)  Gauge  P3 


Time  (msec) 

(b)  Gauge  P2 


Time  (msec) 

(d)  Gauge  P4 


Time  (msec) 


(e)  Gauge  P5 


(f)  Gauge  P6 
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Figure  A12.  (Continued). 


Time  (msec) 

(g)  Gauge  P7 
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(i)  Gauge  P9 
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(h)  Gauge  P8 


G)  Gauge  PIO 


(m)  Gauge  P12 
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Figure  A13.  Comparison  for  CTH  at 
side-on  pressure  gauge  PR. 


BLS  8-ftx8-ft  C2SQ  configuration  with  box 

Calculation  with  SESAME  EOS 

The  CTH  calculation  of  the  8-ftx8-ft  C2SQ  configuration  with  box  used 
473  million  cells  and  required  199  hr  of  wall-clock  time  using  1,024  cores. 
It  required  two  restarts  to  run  the  simulation  out  to  70  msec.  The  primary 
objective  was  to  calculate  the  time-histories  for  the  pressure  gauges 
located  on  the  box  structure  (see  Figure  8(d)).  Since  this  calculation  used 
half  symmetry,  gauges  located  on  the  opposing  half  of  the  box  were 
mirrored  about  the  plane  of  symmetry  when  possible. 

Before  showing  the  time-history  traces,  the  CTH  pressure  fields  between 
20  and  43  msec  are  presented  in  Figure  A14  to  give  the  reader  a  sense  of 
how  the  incident  and  reflected  shocks  transition  with  time.  The  units  of 
pressure  in  the  figure  are  dynes/cm^.  At  20  msec,  the  incident  shock  is 
entering  the  Cascade  section  of  the  BLS  with  a  planar  front  and  a  relatively 
homogeneous  mass  of  air  at  constant  pressure  behind  it  (from 
approximately  8.8  m  to  11.3  m).  A  perturbed  pressure  field  has  engulfed 
the  box  structure  at  25  msec  with  a  reflected  shock  visible  in  front  of  the 
box.  A  planar  reflected  shock  is  visible  moving  upstream  at  31  msec.  At 
37  msec,  the  reflected  shock  front  moving  upstream  has  engulfed  the  box 
structure  and  interacts  with  the  mass  of  air  still  moving  downstream.  This 
interaction  produces  a  nonplanar  front.  A  complex  pressure  field  both 
upstream  and  downstream  of  the  box  structure  is  observed  at  43  msec. 
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Figure  A14.  CTH  pressure  state  plots,  8-ftx8-ft  C2SQ  configuration  with  box. 
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The  CTH  pressure  and  impulse  traces  from  the  box  gauges  are  compared 
to  the  measured  data  in  Figures  A15  to  A18.  The  front  face  gauge  records 
are  plotted  in  Figure  A15,  the  side  gages  in  Figure  A16,  the  top  surface 
gages  in  Figure  A17,  and  the  back  face  gages  in  Figure  A18.  In  addition,  the 
CTH  and  experimental  traces  at  the  P3  and  P9  gage  locations  are  com¬ 
pared  in  Figure  A19.  All  of  the  CTH  time-histories  were  shifted  -4.11  msec 
to  match  the  measured  TOA  at  gauge  PBFi. 


Figure  A15.  Gauge  comparisons  for  CTH,  front  face. 
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Figure  A16.  Gauge  comparisons  for  CTH,  right  face. 
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(b)  Gauge  PBR2 
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Figure  A17.  Gauge  comparisons  for  CTH,  top  face. 
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Figure  A18.  Gauge  comparisons  for  CTH,  back  face. 
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Figure  A19.  Gauge  comparisons  for  CTH,  target  waii. 
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The  ability  of  CTH  to  qualitatively  match  the  measured  data  for  gauges  on 
the  box  varied  from  good  to  poor.  In  general,  very  good  agreement  was 
observed  between  the  CTH  and  the  experimental  traces  up  to  the  arrival  of 
the  second  pressure  peak,  i.e.,  the  arrival  of  the  reflected  shock  from  the 
end  of  the  BLS.  After  that  time  and  for  a  duration  of  approximately 
10  msec,  all  of  the  CTH  traces  from  the  front,  right,  and  top  face  gauge 
locations  fall  below  the  experimental  results.  After  that  time,  the 
magnitude  of  the  pressure  traces  increased  and  exhibited  good  agreement 
with  the  measurements.  The  calculated  pressure  traces  on  the  back  face 
showed  better  agreement  with  the  measured  data  than  the  traces  on  the 
other  surfaces.  Pressure-time-histories  for  two  gauges  on  the  bulkhead  at 
the  end  of  the  BLS  are  plotted  in  Figure  A20.  The  CTH  traces  show  very 
good  agreement  with  the  experimental  results. 

The  incident  pressure  fields  in  the  vicinity  of  the  box  structure  between 
22  and  25  msec  are  shown  in  an  elevation  view  1.0  cm  off  of  the  symmetry 
plane  (Figure  A21)  and  in  a  plan  view  at  the  gauge  PBF3  elevation 
(Figure  A22).  The  top  of  each  image  is  downstream  of  the  box.  The  un- 
shifted  TOA  of  the  incident  shock  on  the  front  face  of  the  box  was 
21.8  msec.  At  22  msec,  the  reflected  shock  from  the  front  face  is  visible  in 
Figure  A2i(a)  and  A22(a).  At  23  msec  (Figures  A2i(b)  and  A22(b)),  the 
reflected  shock  has  progressed  over  the  top  surface  and  the  side  of  the  box 
and  begins  to  wrap  around  the  back  edges. 

The  unshifted  TOA  of  the  incident  shock  at  the  end  of  the  BLS  is 
27.3  msec.  At  31  msec,  the  reflected  shock  is  approximately  halfway 
between  the  back  face  of  the  box  and  the  end  of  the  BLS.  Figures  A23  and 
A24  show  the  reflected  shock  impinging  on  the  bottom,  back-side  of  the 
box  at  33  msec.  By  35  msec,  the  shock  is  wrapping  around  the  top  of  the 
box  onto  the  front  face. 
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Figure  A20.  CTH  pressure  plots  of  incident  shock 
around  box,  elevation  view. 
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Figure  A21.  CTH  pressure  plots  of  incident  shock 
around  box,  plan  view. 
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Figure  A22.  Later-time  CTH  pressure  state  piots. 
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(b)  Pressure  at  35  msec 


Figure  A23.  CTH  pressure  piots  of  reflected  shock 
around  box,  eievation  view. 
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Calculation  with  ideal  gas  EOS 

This  calculation  used  473  million  cells,  required  180  hr  of  wall-clock  time 
using  1,024  cores,  and  required  two  restarts  to  run  the  simulation  out  to 
70  msec.  The  calculation  was  conducted  to  see  if  the  ideal  gas  EOS 
provided  improved  pressure  and  impulse  traces  from  the  box  gauges.  The 
time-histories  were  shifted  -2.91  msec  to  match  the  TOA  at  gauge  PBFi. 
Figure  A24  compares  the  traces  from  three  gauge  locations,  one  on  each  of 
the  front,  side,  and  back  of  the  box  structure.  On  the  left  side  of  the  figure, 
the  CTH  waveforms  are  compared  to  the  experimental  results  for  gauges 
PBFi,  PBRi,  and  PBBi,  while  only  the  CTH  waveforms  are  shown  on  the 
right  side. 

At  the  PBFi  gauge  location  (Figures  A24  (a)  and  A24  (b)),  the  CTH 
impulse  traces  were  nominally  the  same  up  to  36  msec.  Between  36  and 
45  msec,  the  ideal  gas  impulse  trace  plots  above  the  SESAME  trace  and 
after  45  msec  it  lies  below.  The  simulated  impulse  traces  were  well  below 
the  experimental  results.  The  simulated  pressure  and  impulse  records  at 
the  PBRi  gauge  location  (Figures  A24  (c)  and  A24(d))  matched  the 
experimental  results  well  until  37  msec.  After  that  time,  the  SESAME 
pressures  drop  well  below  the  experimental  and  a  similar  drop  in  pressure 
is  observed  in  the  ideal  gas  pressures  at  42  msec.  Both  simulation  pressure 
traces  track  each  other  after  60  msec.  At  the  PBRi  gauge  location 
(Figures  A24(e)  and  A24(f)),  the  simulated  pressure  and  impulse  records 
track  the  experimental  results  well  until  40  msec  at  which  time  the 
SESAME  traces  fall  below  both  the  ideal  gas  and  experimental  traces.  The 
ideal  gas  impulse  record  begins  to  diverge  from  the  experimental  at 
50  msec. 

Overall,  the  ideal  gas  EOS  provided  a  marginal  improvement  over  the 
SESAME  EOS.  Some  records  exhibited  an  improved  match  to  the  experi¬ 
mental  data,  and  others  showed  a  degraded  comparison. 


Overpressure  (psi)  Overpressure  {psi)  Overpressure  (psi) 


ERDC/GSL  TR-16-31 


77 


Figure  A24.  Gauge  comparisons  between 
SESAME  and  ideai  gas  EOS. 


(a)  Gauge  PBFl  (b)  Gauge  PBFl;  CTH  only 


(c)  Gauge  PBRl  (d)  Gauge  PBRl;  CTH  only 


(e)  Gauge  PBBl 


(f)  Gauge  PBBl:  CTH  only 
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Lessons  learned 

The  steel  walls  of  the  BLS  were  originally  modeled  in  CTH  with  a  rigid 
material.  This  should  have  maximized  the  computational  time-steps  since 
wave  speeds  in  air  would  be  much  larger  than  steel  or  any  other  metal. 
Early  3-D  calculations  with  the  rigid  walls  suffered  from  an  exponentially 
decaying  time-step.  For  a  5-msec  run,  the  time-step  decayed  from  1.0  to 
0.2  psec.  Simulations  of  just  the  driver  and  the  CVC  sections  of  the  BLS 
were  taking  48  hr  of  wall-clock  time  using  1,024  cores  for  an  11  msec  run. 
Various  CTH  input  options  were  tried  to  mitigate  the  decaying  time-steps 
without  success,  i.e.,  changing  from  SESAME  to  the  ideal  gas  EOS,  using 
different  multiple  material  and  temperature  models,  and  removing  the 
striker  from  the  geometry.  The  cause  of  the  problem  was  finally  recognized 
after  temperatures  on  the  order  of  3,100  deg  Kelvin  (0.27  eV)  were 
observed  in  the  region  of  the  striker  (Figure  A25).  The  problem  was  solved 
by  removing  the  rigid  material  and  replacing  it  with  a  non-rigid  material. 
Given  the  1.27-cm-thick  walls  and  the  nominal  0.5-cm-sized  cells  in  the 
mesh,  the  number  of  cells  through  the  wall  thickness  was  insufficient  for  a 
rigid  material.  Unfortunately,  this  problem  severely  impacted  the  progress 
of  the  CTH  modeling  effort  as  it  required  a  significant  amount  of  time  to 
correct. 


Figure  A25.  Temperatures  in  the  region  of  the  striker. 
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Summary 

CTH  was  used  to  simulate  the  blast  environment  for  three  configurations 
of  the  ERDC  BLS.  Qualitatively,  CTH  was  able  to  capture  the  system 
responses  for  each  configuration.  The  calculated  time-histories  for  both 
the  GSA  and  the  8 -ft  x  8 -ft  C2SQ  configuration  matched  the  experimental 
records  very  well.  For  the  8-ftx8-ft  C2SQ  configuration  with  box  structure, 
CTH  underpredicted  the  impulse  for  all  gauges  on  the  box. 

The  combination  of  a  large  length  scale  and  thin  BLS  walls  did  not  fit  well 
into  the  strengths  of  CTH.  It  is  felt  the  code  would  perform  better  for  prob¬ 
lems  involving  single  or  multiple  structures  in  an  open-air  environment 
subjected  to  blast  from  an  explosive  source.  CTH  has  a  unique  capability  to 
handle  material  strength  in  its  Eulerian  solution  approach.  This  capability 
is  not  available  in  typical  CFD  codes  and  would  be  necessary  if  any  simu¬ 
lated  structures  exhibited  elastic  or  inelastic  deformations. 
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Appendix  B;  DYSMAS  Modeling 

DYSMAS  background 

The  coupled  code  DYSMAS  (Harris  et  al.  2014)  was  developed  by  NSWC/ 
Indian  Head  to  simulate  underwater  explosion  (UNDEX)  environments 
and  the  ensuing  marine  vessel  response  resulting  from  the  fluid-structure 
interaction.  The  loosely  coupled  solution  approach  has  three  components: 
an  Eulerian  code  Gemini  (Wardlaw  et  al.  2003)  that  performs  the  fluid 
flow  calculation,  a  Lagrangian  code  ParaDyn  (DeGroot  et  al.  2013)  that 
performs  the  structural  response  calculation,  and  a  coupler  that  transfers 
information  between  Gemini  and  ParaDyn.  This  enables  DYSMAS  to  exe¬ 
cute  calculations  that  include  explicit  shock  fronts,  bubble  jets,  structural 
failure,  and  fluid  breakthrough.  Included  in  Gemini  are  the  capabilities  to 
handle  compressible  and  incompressible  fluids.  Gemini  employs  a  time- 
split,  second-order  Godunov  scheme  to  solve  the  Euler  equations.  ParaDyn 
is  a  3-D  explicit  finite  element  program  for  analyzing  the  dynamic 
response  of  solids  and  structures. 

DYSMAS  is  a  fluid-structure  interaction  set  of  codes  that  was  designed 
especially  for  weapon  lethality  studies,  where  the  explosive  charge  is 
detonated  in  close  proximity  to  the  target  structure.  For  such  cases,  there 
is  strong  interaction  between  the  UNDEX  loading  and  the  target  structure, 
where  the  motion  of  the  structure  significantly  modifies  the  loading.  The 
suite  of  codes  can  be  applied  to  any  fluid-structure  interaction  problem  or 
the  individual  portions  used  to  simulate  their  regime  alone  (e.g.,  Gemini 
applied  to  a  fluid  problem).  DYSMAS  has  been  under  continuous  develop¬ 
ment  as  part  of  a  U.S.-German  Project  Agreement  entitled  “Enhanced 
Undersea  Weapons  Effectiveness  and  Ship  Survivability  through  the 
Application  of  Validated  Computer  Codes.” 

The  Gemini  code  is  an  Eulerian  solver  for  the  fluid  equations  of  motion 
and  was  designed  specifically  to  simulate  explosions  in  water;  however,  it 
has  been  applied  to  other  environments,  including  explosions  in  air  and 
soil.  The  ParaDyn  code  is  an  explicit  dynamics  Lagrangian  finite  element 
code  that  solves  the  structural  equations  of  motion.  The  existing  interface 
routines  allow  Gemini  and  ParaDyn  to  exchange  information  in  order  for 
the  fluid  and  structural  integration  to  advance  in  time.  This  Eulerian- 
Lagrangian  combination  provides  a  powerful  capability  that  has  been 
extensively  validated  for  complex  UNDEX  phenomena  including  explosive 
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shock  propagation,  bubble  formation  and  jetting,  and  fluid-structure 
interaction.  Gemini  uses  a  Cartesian  structured  mesh  (i.e.,  fixed  rectangu¬ 
lar  cells)  and  supports  setups  in  varied  dimensions  (i.e.,  i-D,  2-D,  or  3-D). 

ParaDyn  is  the  parallel  implementation  of  the  finite  element  code  DYNA3- 
D  (Whirley  and  Englemann  1993).  DYNA3-D  is  based  on  a  finite  element 
discretization  of  the  three  spatial  dimensions  and  a  finite  difference 
discretization  of  time.  The  explicit  central  difference  method  is  used  to 
integrate  the  equations  of  motion  in  time.  The  central  difference  method  is 
conditionally  stable,  and  stability  is  governed  by  the  Courant  limit  on  the 
time-step.  DYNA3-D  uses  a  lumped  mass  formulation  for  efficiency.  This 
produces  a  diagonal  mass  matrix,  which  renders  the  solution  of  the 
momentum  equation  trivial  at  each  step  in  that  no  simultaneous  system  of 
equations  must  be  solved.  The  contribution  of  the  element  reaction  forces 
and  applied  loads  are  accumulated  into  a  nodal  force  vector,  which  is  then 
divided  by  the  lumped  nodal  mass  to  get  the  nodal  acceleration.  This  is 
then  integrated  forward  in  time  using  a  central  difference  approximation 
to  get  the  nodal  velocity  and  displacement. 

The  ParaDyn  code  passes  a  list  of  “interface  surfaces”  to  Gemini  at  the 
start  of  a  coupled  calculation  that  describes  the  location  of  the  structure 
within  the  fluid  domain  (i.e.,  the  wetted  surface).  Each  discrete  surface  is 
defined  solely  by  three  or  four  nodes  specified  in  the  “interface  element” 
section  of  the  ParaDyn  input  deck.  Loads  are  passed  from  Gemini  to 
ParaDyn,  inducing  a  response  and  potential  motion  of  the  structure.  That 
new  location  and  state  of  the  structure  are  then  passed  back  to  Gemini  for 
further  loading.  This  coupled  treatment  occurs  at  every  time-step  in  the 
calculation. 

Modeling  overview 

All  of  the  DYSMAS  calculations  were  run  on  the  Air  Force  Research 
Laboratory  (AFRL)  HPC  platform  Spirit.  Spirit  has  4,590  compute  nodes 
with  16  cores  per  nodes.  Each  compute  node  has  a  core  speed  of  2.6  GHz 
and  32  GB  of  accessible  memory. 

Full  3-D  calculations  were  conducted  for  all  cases  considered  in  the  BLS 
modeling  effort.  The  BLS  structure  was  modeled  as  a  Lagrangian  con¬ 
struct,  which  was  then  inserted  into  the  Eulerian  air  domain.  The  BLS 
structure  was  composed  of  shell  elements  with  a  nominal  thickness  of 
0.5  in.  The  nodal  displacements  were  fixed  so  that  the  structure  would  act 
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in  a  completely  non-responding  fashion.  A  typical  model  is  shown  in 
Figure  Bi.  Transmitting  boundary  conditions  were  applied  to  all  mesh 
boundaries  in  the  Gemini  portion  of  the  calculation.  The  Eulerian  domain 
extended  a  moderate  distance  away  from  the  BLS  structure.  A  uniform  cell 
size  was  specified  throughout  the  Eulerian  domain  encompassed  by  the 
BLS.  A  graded  mesh  was  specified  beyond  the  extents  of  the  BLS.  The 
driver  gas  was  represented  by  a  volume  of  high  pressure  air  that  began  to 
propagate  a  time  zero  (see  Figure  B2).  The  diaphragm  and  its  failure  were 
not  modeled.  An  ideal  gas  EOS  was  used  to  model  both  the  ambient  air 
and  the  high  pressure  air  contained  within  the  driver  section.  These  were 
modeled  as  two  separate  Eulerian  materials  in  the  problem  setup.  The 
ambient  air  was  assigned  an  initial  pressure  and  temperature  of  14.5  psia 
and  77°F,  respectively.  The  driver  air  was  assigned  an  initial  pressure  of 
1,298  psig  and  initial  temperature  of  77°F.  The  BLS  structure  was  numeri¬ 
cally  rigid  and  fixed  in  space  so  as  to  not  respond/move.  When  appropri¬ 
ate,  the  striker  and  grill  were  explicitly  modeled  within  ParaDyn. 

A  series  of  calculations  were  performed  to  assess  the  effects  of  the  striker 
and  grill  on  the  flow  environment.  A  simulation  was  performed  without 
either  of  those  constructs,  one  with  just  the  striker  and  one  with  both  the 
striker  and  the  grill.  The  striker  and  grill  were  modeled  with  solid  elements 
within  ParaDyn  and  a  coupling  interface  applied  to  their  surfaces  to  interact 
with  the  fluid  flow.  Both  are  depicted  in  Figure  B3,  shown  without  the  BLS 
structure.  Figure  B4  shows  the  flow  field  for  all  three  configurations.  The 
top  row  depicts  the  configuration  without  the  striker  or  grill;  the  second  row 
of  images  depicts  the  striker  only;  and  the  images  on  the  third  row  include 
both  the  striker  and  grill.  Qualitatively,  there  is  little  difference  in  the 
pressure  field,  suggesting  that  the  striker  and  grill  have  only  a  minor  effect 
on  the  shock  and  trailing  high  density  gas  that  propagate  down  the  BLS.  A 
quantitative  comparison  is  provided  in  Figure  B5,  which  depicts  the  pres¬ 
sure  history  calculated  near  the  center  of  the  target  wall.  Although  each 
construct  has  an  effect  on  the  waveform,  the  effects  are  still  minor.  An 
approximate  5  percent  reduction  in  peak  pressure  is  seen  from  each  part, 
and  at  30  msec,  the  striker  reduces  the  impulse  by  about  5  percent  and  the 
grill  an  additional  2  percent.  The  arrival  times  are  delayed  by  each  part, 
approximately  0.2  msec  each  for  both  the  striker  and  grill. 
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Figure  Bl.  Typical  DYSMAS  model  (8-ftx8-ft 
C2SQ  configuration  with  box  shown). 
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Figure  B3.  Depiction  of  striker  and  griii  modeis. 


Figure  B4.  Effect  of  striker  and  griii  modeis  on  the  flow  fieid. 


Figure  B5.  Effects  of  striker  and  griii. 
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Another  set  of  calculations  was  performed  on  one  of  the  initial  models 
with  an  increasing  mesh  density  within  the  fluid  domain  to  assess  mesh 
convergence.  If  grid  convergence  is  obtained,  then  solutions  from  the  two 
finest  meshes  should  be  identical,  and  there  is  confidence  that  no  further 
refinement  is  needed.  Three  mesh  densities  were  used,  having  a  minimum 
cell  size  of  2.0, 1.0,  and  0.5  cm,  respectively.  This  led  to  fluid  meshes  of 
approximately  18.3  million,  135  million,  and  1.0  billion  cells.  The  struc¬ 
tural  mesh  was  unchanged  throughout  this  exercise,  as  altering  its  mesh 
density  should  have  very  minimal  effects  on  the  flow  field  and  behavior 
due  to  it  being  rigid  and  fixed  in  space. 

Run  time  statistics  are  provided  in  Table  Bi.  Results  from  the  convergence 
study  are  provided  in  Figure  B6,  which  depicts  the  pressure  history  and 
integrated  impulse  at  a  single  location  on  the  target  wall  (gauge  P5).  These 
results  are  typical  of  those  at  the  other  gauge  stations.  Refining  the  mesh 
to  a  0.5-cm  cell  size  produces  a  slightly  different  waveform  when  observed 
locally,  but  the  global  waveforms  match  well.  There  is  a  0.85  percent 
difference  in  total  impulse  at  125  msec,  indicating  the  coarser  mesh  with  a 
1.0-cm  cell  size  is  sufficient  to  resolve  the  pressures  imparted  to  the 
structure  within  this  configuration. 


Table  Bl.  DYSMAS  run  time  statistics  from  convergence  study. 


Mesh 

Resolution 

(cm) 

Analysis  End 
Time  (msec) 

No.  Gemini 
Time-steps 

No. 

Processors 

CPU  Time 
(hours) 

Approximate 

Processor- 

hours 

2 

150.0 

10,252 

512 

8.51 

4,357 

1 

125.7 

14,418 

1024 

16 

16,384 

0.5 

150.0 

28,335 

4048 

44.13 

178,638 
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Figure  B6.  Convergence  study  results. 


(a)  Global  comparisons 


(b)  Early  time 


GSA  configuration,  caiibration  tests  30a  and  30b 

The  DYSMAS  model  of  the  GSA  configuration  is  depicted  in  Figure  By.  All 
of  the  BLS  is  explicitly  modeled,  including  the  vents,  Cascade,  and  target 
wall.  This  test  included  the  striker  and  grill,  and  this  configuration  is 
captured  within  the  model  (see  Figure  B8). 

Figure  B7.  GSA  configuration  model  with  slice 
displaying  initial  pressure  state. 
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Figure  B8.  GSA  configuration  modei,  ciose-up 
of  striker  and  driver  pressure  state. 
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A  constant  cell  size  of  i.o  cm  was  used  throughout  the  interior  of  the  BLS 
model.  On  the  exterior,  a  graded  mesh  extended  to  ±50.0  cm  in  the  x  and 
y  directions  and  ±20.0  cm  in  the  z  direction  along  the  length  of  the  BLS 
model.  This  resulted  in  a  mesh  of  approximately  136  million  cells.  A  cross 
section  of  a  portion  of  the  grid  superimposed  on  a  contour  plot  of  the 
initial  pressure  is  provided  in  Figure  B9.  The  structure  model  consists  of 
29,128  shell  elements  in  the  superstructure  of  the  BLS  and  3,886  solid 
elements  within  the  striker  and  grill. 


Figure  B9.  GSA  configuration  modei  mesh. 
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A  series  of  pressure  state  plots  depicting  the  shock  and  trailing  high 
density  gas  flow  is  provided  in  Figure  Bio.  They  depict  the  gas  leaving  the 
driver,  propagating  down  the  BLS,  venting  from  the  CVC  section,  reflect¬ 
ing  off  the  target  wall  at  the  downstream  end,  and  flowing  back  upstream. 
The  time  stamps  in  the  individual  state  plots  are  given  in  units  of  seconds. 
The  shock  propagating  laterally  from  the  driver  section  at  1.9  msec  is 
caused  by  the  interaction  of  the  high  pressure  gas  with  the  surrounding 
Lagrangian  structure  at  startup  of  the  calculation.  While  the  initial  shock 
propagates  in  a  very  planar  fashion,  there  is  a  considerable  amount  of 
heterogeneity  in  the  flow  field  behind  the  front,  especially  at  the  later 
times  (e.g.,  see  state  plots  at  27.5  and  40.0  msec). 

A  comparison  of  the  calculated  pressure  history  at  select  gauges  located  on 
the  target  wall  is  provided  in  Figure  Bn.  The  results  for  gauges  P6,  P9,  and 
Pro  are  not  included  since  they  are  at  repeat  locations  from  the  center 
point.  All  of  the  gauges  exhibit  the  same  basic  pattern  with  small  oscilla¬ 
tions  around  the  global  shock  pattern.  The  most  obvious  differences  occur 
early-on  in  the  reflected  shock  history;  however,  these  differences  are 
relatively  small  indicating  that  DYSMAS  is  predicting  a  planar  front 
interacting  with  the  target  wall.  Gauge  P3,  which  resides  at  the  center  of 
the  target  wall,  exhibits  the  greatest  fluctuation.  The  fluctuations  in  the 
pressure  history  generally  get  smoother  as  one  moves  away  from  the 
centerline. 

Comparisons  between  the  measured  and  calculated  data  are  provided  in 
Figure  B12  for  several  gauge  locations.  The  calculation  was  run  to  63  msec, 
which  was  sufficient  to  capture  the  first  pulse  in  the  measured  data,  as  well 
as  the  initial  rise  for  the  second.  The  calculated  results  compare  well  with 
the  measured  data  at  all  gauge  locations.  Note,  the  calculated  results  are 
“as-is”  and  have  not  been  shifted  in  time  to  match  the  measured  shock 
arrival. 
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Figure  BIO.  DYSMAS  pressure  state  plots  for  GSA  configuration. 


Figure  Bll.  Calculated  gauge  pressure,  GSA  configuration. 
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Figure  B12.  Gauge  comparisons  for  DYSMAS,  GSA  configuration. 
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Figure  B12.  (Continued). 
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A  i.o-cm  uniform  cell  size  was  specified  for  the  region  within  the  BLS,  with 
a  graded  Eulerian  mesh  utilized  outside.  The  analysis  end  time  was 
62.7  msec.  It  took  16  hr  of  CPU  time  and  1,016  processors  to  run  this 
calculation  to  62.7  msec.  In  the  calculation,  8  processors  were  alloted  to 
ParaDyn  and  1,008  to  Gemini.  There  were  approximately  136  million  cells 
in  the  Gemini  mesh.  It  took  13,370  time-steps  in  the  Gemini  portion  of  the 
coupled  calculation  to  run  the  analysis  to  62.7  msec. 

Calibration  tests  with  8-ftx8-ft  C2SQ  configuration 

The  DYSMAS  model  for  the  second  case  included  the  striker,  grill,  and 
extended  section  of  the  BLS.  The  Eulerian  domain  contained 
approximately  165.5  million  cells,  while  the  BLS  structure  consisted  of 
29,128  shells  and  3,886  solid  elements.  This  configuration  and  a  cross 
section  of  the  pressure  field  at  an  early  time  are  shown  in  Figure  B13.  A 
progression  through  time  of  the  pressure  state  is  provided  in  Figure  B14, 
which  encompasses  the  high-density  gas  exiting  the  driver  early-on,  the 
shock  and  trailing  high-density  gas  propagating  down  the  BLS,  reflection 
off  of  the  target  wall,  and  propagation  back  upstream  toward  the  vents. 

The  time  stamp  on  each  figure  is  given  in  units  of  seconds. 

Similar  to  the  first  case,  this  calculation  exhibits  a  clean,  planar  shock 
front  and  more  chaotic  flow  behavior  behind  the  front.  By  observing  the 
particle  velocity  along  the  length  of  BLS  and  the  density,  this  heterogene¬ 
ous  behavior  can  be  clearly  seen  (Figure  B15).  The  contact  surface,  which 
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Figure  B13.  Early  time  pressure  state,  8-ftx8-ft 
C2SQ  configuration. 


Figure  B14.  Pressure  state  plots  for  8-ftx8-ft  C2SQ  configuration. 
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Figure  B15.  Velocity  and  density  state  plots  at  80  msec. 


(a)  Flow  velocity 


(b)  Density 


is  the  interface  between  the  driver  and  downstream  gas,  follows  the  shock 
front  but  only  propagates  a  limited  distance  down  the  tube.  The  figure 
indicates  that  the  contact  surface  extends  downstream  to  the  end  of  the 
Cascade,  but  does  not  go  further  (determined  from  additional  state  plots, 
not  shown);  however,  the  shock  front  continues  onwards  as  a  disturbance 
propagating  through  the  remaining  ambient  air. 

The  calculated  pressure  histories  at  various  locations  on  the  target  wall  are 
shown  in  Figure  Bi6.  The  mostly  planar  nature  of  the  shock  front  is 
present,  with  only  slight  variations  noted  in  the  early  reflected  shock 
illustrated  in  Figure  Bi6(b).  Comparisons  between  the  measured  and 
calculated  data  for  the  target  wall  gauges  are  provided  in  Figure  B17.  The 
calculated  results  are  presented  “as-is”  and  have  not  been  shifted  in  time 
to  match  the  measured  shock  arrival.  Overall,  one  notes  good  correlation 
between  the  calculated  and  measured  waveform,  with  a  slight  tendency  to 
overpredict  the  arrival  time  (about  1.3  msec  later).  The  only  significant 
disparities  are  for  gauge  P7;  however,  the  results  from  this  gauge  are 
suspect  when  compared  with  the  amplitude  and  general  structure  of  the 
waveform  at  the  other  gauge  locations. 
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Figure  B16.  Calculated  gauge  pressure,  8-ftx8-ft 
C2SQ  configuration. 


(a)  Global  comparisons  (b)  Early  time 


Figure  B17.  Gauge  comparisons  for  DYSMAS,  8-ftx8-ft 
C2SQ  configuration. 
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Figure  B17.  (Continued). 
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(f)  Gauge  P6 
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(h)  Gauge  P8 
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Figure  B17.  (Continued). 


(k)  Gauge  Pll  (i)  Gauge  P12 


Comparisons  for  the  side-on  gauge  PR  are  provided  in  Figure  Bi8.  The 
DYSMAS  results  have  been  time-shifted  by  -2.44  msec  to  match  the 
measured  arrival  time.  This  was  done  to  ease  comparison  of  the  measured 
and  calculated  waveforms.  One  notes  that  very  good  correlation  is  noted 
between  the  two. 
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Figure  B18.  Comparisons  for  gauge 
PR,  8-ftx8-ft  C2SQ  configuration. 


A  i.o-cm  uniform  cell  size  was  specified  for  the  region  within  the  BLS,  with 
a  graded  Eulerian  mesh  utilized  outside.  The  analysis  end  time  was 
150  msec.  It  took  31.85  hr  of  CPU  time  and  1,016  processors  to  run  this 
calculation  to  150  msec.  In  the  calculation,  8  processors  were  alloted  to 
ParaDyn  and  1,008  to  Gemini.  There  were  approximately  136  million  cells 
in  the  Gemini  mesh.  It  took  19,212  time-steps  in  the  Gemini  portion  of  the 
coupled  calculation  to  run  the  analysis  to  150  msec. 

BLS  8-ftx8-ft  C2SQ  configuration  with  box 

This  BLS  configuration  is  the  same  as  in  the  previous  case  with  the  excep¬ 
tion  of  a  box  structure  inserted  into  the  C2SQ  section  (see  Figure  B19).  As 
with  the  previous  case,  the  BLS  model  included  both  the  striker  and  grill 
(see  Figure  B20).  An  additional  view  looking  upstream  from  the  calibra¬ 
tion  plate  towards  the  driver  section  is  provided  in  Figure  B21. 

A  time  sequence  of  the  pressure  propagation  is  shown  for  the  pressure 
wrapping  around  the  structure  in  Figures  B22  and  B23,  picturing  the  pres¬ 
sure  inside  and  overhead  views,  respectively.  The  velocity  and  density 
plots  in  Figure  B24  illustrate  that  the  contact  surface  (interface  between 
the  driver  gas  and  downstream  air)  only  extends  to  the  end  of  the  Cascade. 
Thus,  the  pressure  history  for  gauges  residing  on  the  box  should  be  fairly 
clean  and  more  consistent,  as  opposed  to  the  box  being  placed  at  a  location 
further  upstream.  The  lack  of  noise  in  the  measured  data  bears  this  out. 
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Figure  B19.  DYSMAS  model  of  8-ftx8-ft  C2SQ  configuration  with  box. 


Figure  B20.  Close-up  of  striker  and 
grill  at  the  driver  end. 


Figure  B21.  End-on  view  looking  upstream 
to  driver  section. 
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Figure  B22.  Pressure  state  plots,  8-ftx8-ft  C2SQ 
configuration  with  box,  elevation  view. 


The  calculated  pressure  histories  for  gauges  on  the  box  structure  are  out¬ 
lined  in  Figures  B25  through  B29.  These  figures  provide  a  global  compari¬ 
son  of  the  calculated  results  along  with  a  comparison  of  the  early  time 
shock  interaction  for  data  on  each  of  the  faces  (except  the  right  face,  which 
is  comparable  to  that  of  the  left  one).  In  general,  one  notes  good  consis¬ 
tency  in  the  calculated  results  for  gauges  on  each  face.  There  are  slight 
differences  for  the  left  and  right  face  gauges  due  to  grid  alignment;  how¬ 
ever,  these  are  negligible.  The  calculated  pressures  on  each  face  provide 
insight  into  the  shock  interaction  with  the  box  structure.  The  initial  shock 
distinctly  strikes  the  front  face  of  the  box  and  continues  propagating  for¬ 
ward  in  time,  wrapping  around  the  structure  as  seen  in  the  state  plots 
above.  Once  the  shock  reflects  off  the  downstream  target  wall,  it  again 
wraps  around  the  structure  with  a  large  reflected  pressure  on  the  back  sur¬ 
face  and  refracted  pressures  on  the  remaining  surfaces.  This  process  is 
also  shown  in  Figure  B30,  where  pressures  for  the  front,  back,  left  side. 
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and  top  are  plotted  against  each  other.  The  relative  timings  are  correct  for 
both  the  initial  and  reflected  shocks,  as  are  the  relative  peaks  with  the 
highest  being  the  first  face  encountered  (front  face  for  the  initial  shock  and 
back  face  for  the  reflected  shock)  after  which  the  refracted  shock  front  has 
lower  peaks  on  the  side  and  back  faces. 


Figure  B23.  Pressure  state  plots,  8-ftx8-ft  C2SQ 
configuration  with  box,  plan  view. 
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Figure  B24.  Velocity  and  density  state  plots  at  50  msec. 


(a)  Flow  velocity 


(b)  Density 


Figure  B25.  Calculated  gauge  pressure  for  DYSMAS,  front  face. 


(a)  Global  comparisons 


(b)  Early  time 
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Figure  B26.  Calculated  gauge  pressure  for  DYSMAS,  back  face. 


(a)  Global  comparisons  (b)  Early  time 


Figure  B27.  Calculated  gauge  pressure  for  DYSMAS,  left  face. 


(a)  Global  comparisons 


(b)  Early  time 
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Figure  B28.  Calculated  gauge  pressure  for  DYSMAS,  right  face. 


(a)  Global  comparisons  (b)  Early  time 


Figure  B29.  Calculated  gauge  pressure  for  DYSMAS,  top  face. 
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(b)  Early  time 
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Figure  B30.  Comparison  of  caicuiated  gauge  pressure  on  each  face. 


(a)  Giobai  comparisons 


(b)  Eariytime 


Comparisons  between  the  measured  and  calculated  data  for  the  gauges  on 
the  box  structure  are  provided  in  Figures  B31  through  B35.  The  calculated 
results  are  presented  “as-is”  and  have  not  been  shifted  in  time  to  match 
the  measured  shock  arrival.  Overall,  one  notes  good  correlation  between 
the  calculated  and  measured  waveforms  with  a  slight  tendency  to  over¬ 
predict  the  arrival  time  (about  1.8  msec  later).  Similar  comparisons  are 
provided  for  the  target  wall  (see  Figure  B36).  Again,  good  correlation  is 
noted  between  the  calculated  and  measured  data. 

A  1.0-cm  uniform  cell  size  was  specified  for  the  region  within  the  BLS, 
with  a  graded  Eulerian  mesh  utilized  outside.  The  analysis  end  time  was 
150  msec.  It  took  37.48  hr  of  CPU  time  and  1,016  processors  to  run  this 
calculation  to  150  msec.  In  the  calculation,  8  processors  were  alloted  to 
ParaDyn  and  1,008  to  Gemini.  There  were  approximately  165.5  million 
cells  in  the  Gemini  mesh.  It  took  19,600  time-steps  in  the  Gemini  portion 
of  the  coupled  calculation  to  run  the  analysis  to  150  msec. 


Overpressure  (psi)  Overpressure  (psi)  Overpressure  (psi) 


ERDC/GSL  TR-16-31 


105 


Figure  B31.  Gauge  comparisons  for  DYSMAS,  front  face. 
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Figure  B32.  Gauge  comparisons  for  DYSMAS,  back  face. 
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Figure  B33.  Gauge  comparisons  for  DYSMAS,  ieft  face. 


Time  (msec) 


(a)  Gauge  PBLl 


(b)  Gauge  PBL2 


(c)  Gauge  PBL3 


(d)  Gauge  PBL4 


Impulse  (psi-msec)  Impulse  (psi-msec) 


Overpressure  (psi)  Overpressure  (psi) 


ERDC/GSL  TR-16-31 


108 


Figure  B34.  Gauge  comparisons  for  DYSMAS,  right  face. 
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Figure  B35.  Gauge  comparisons  for  DYSMAS,  top  face. 
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Figure  B36.  Gauge  comparisons  for  DYSMAS,  target  waii. 
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Appendix  C:  LOCI-Blast  Modeling 

Loci/BLAST  Background 

Loci/BLAST  is  a  modern  computational  fluid  dynamics  (CFD)  code 
designed  to  provide  accurate  simulations  for  a  variety  of  blast  scenarios. 
Loci/BLAST  is  based  on  the  Loci/CHEM  chemically  reacting  flow  code 
developed  at  Mississippi  State  University  (MSU)  by  Dr.  Edward  Luke  and 
his  co-workers  (Luke  1999;  Luke  and  George  2005;  Luke  and  Cinella 
2007).  Loci/CHEM  has  been  in  constant  development  for  more  than 
15  years  and  is  widely  used  at  the  National  Aeronautics  and  Space 
Administration  (NASA)  and  in  the  aerospace  industry  for  CFD  analysis  of 
aircraft  and  launch  vehicle  aerodynamics  and  chemically  reacting  flows  in 
rocket  motors.  Development  of  Loci/BLAST  was  initiated  under  the 
U.S.  Army  Tank  Automotive  Research  and  Development  Engineering 
Center  (TARDEC)  Simulation  Based  Reliability  and  Safety  (SimBRS) 
program  to  provide  a  modern  computational  tool  for  analysis  of 
underbody  blast  effects  on  military  vehicles  (Thompson  et  al.  2010, 

2012a).  Loci/BLAST  was  extended  to  perform  simulations  of  urban  blast 
problems  under  Department  of  Homeland  Security  sponsored  Southeast 
Region  Research  Initiative  (SERRI)  funding  (Thompson  et  al.  2012b). 

Both  Loci/CHEM  and  Loci/BLAST  have  undergone  extensive  verification 
and  validation  testing  over  the  years  and  have  a  proven  history  of 
providing  accurate  and  cost-effective  solutions  for  a  variety  of  CFD  and 
blast  related  problems. 

Loci/BLAST  has  a  wide  range  of  capabilities  that  support  large-scale  blast 
analyses.  Loci/BLAST  is  an  unstructured  cell-centered  finite  volume  code 
that  can  utilize  generalized  grids  consisting  of  a  variety  of  polyhedral  cell 
types  (tetrahedra,  hexahedra,  prisms,  etc.).  Unlike  some  other  codes  used 
in  this  research,  Loci/BLAST  is  fully  Eulerian  (no  operator  splitting)  and 
uses  well-defined,  modern  unstructured  CFD  analysis  methods  such  as  the 
Harten- Lax-van  Leer-Enfeldt  (HLLE)  approximate  Riemann  solver  to 
compute  fluxes  at  cell  faces  and  least-squares  reconstruction  coupled  with 
a  nodal  Baldwin-Barth  limiter  to  provide  second-order  spatial  accuracy. 
Explicit  time  integration  is  performed  using  a  two-stage  second-order 
Total  Variation  Diminishing  (TVD)  Runge-Kutta  scheme.  The  Loci  frame¬ 
work  on  which  both  Loci/BLAST  and  Loci/CHEM  are  based  provides  facil¬ 
ities  for  automatic  parallelization  and  cache  memory  use  optimization  that 
makes  codes  developed  using  the  framework  highly  scalable  and 
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computationally  efficient  in  parallel  computing  environments.  Loci/ 
CHEM  has  demonstrated  scalability  to  10,000  processor  class  jobs  and 
will  scale  to  larger  processor  counts  with  little  or  no  additional  code  modi¬ 
fications.  Other  capabilities  of  Loci/BLAST  that  support  large-scale  blast 
analysis  include  varied  EOS  options  for  energetic  materials,  a  multi-phase 
soil  model,  particulate  transport  models,  and  both  prescribed  initiation 
and  reactive  burn  detonation  models.  Loci/BLAST  can  also  utilize  a 
unique  overset  grid  capability  to  reduce  the  work  load  associated  with 
generating  global  grid  systems  that  contain  a  variety  of  solid  objects  of 
varying  shapes  and  dimensions  such  as  those  found  in  urban  blast 
scenarios. 

The  Loci  framework  provides  a  unique  code  development  paradigm  based 
on  principles  taken  from  relational  database  design  and  artificial  intelli¬ 
gence  instead  of  traditional  procedural  or  object-oriented  approaches. 
Various  analysis  tasks  can  be  instantiated  as  a  set  of  rules  that  operate  on 
prescribed  data  structures  to  generate  a  desired  outcome.  Codes  can  then 
be  rapidly  developed  by  compositing  rules  into  a  set  of  queries  that  per¬ 
form  the  desired  tasks  in  the  correct  order.  This  reduces  the  time  required 
to  generate  complex  multi-physics  codes  because  the  framework  will 
detect  any  discrepancies  in  data  transfers  between  code  components  and 
schedule  computational  tasks  to  execute  in  the  most  computationally 
efficient  manner.  Therefore,  the  time  required  to  implement  future  code 
modifications  required  to  support  the  next  phases  of  the  current  research 
will  in  most  cases  be  shorter  than  with  traditional  programming 
approaches. 

Analysis  overview 

All  of  the  calculations  of  the  ERDC  BLS  tests  using  Loci/BLAST  were  per¬ 
formed  by  Clayton  Mord,  a  graduate  student  in  Aerospace  Engineering  at 
MSU  under  direction  of  Dr.  Richard  Weed  and  constitutes  the  bulk  of  his 
master’s  thesis  research  (Mord  2015).  As  with  the  analyses  performed  with 
the  other  codes  used  in  this  research,  an  initial  set  of  2-D  axisymmetric 
runs  were  made  with  the  initial  BLS  configuration  to  guide  problem  setup 
and  verify  the  code  could  produce  a  planar  shock  wave  in  the  test  section. 
These  runs  were  made  with  a  clean  flow  path  that  did  not  include  the  BLS 
striker  mechanism  and  the  catch  grill.  Subsequent  calculations  were  run 
with  a  variety  of  mesh  configurations  until  the  final  grid  systems  used  for 
the  results  presented  in  the  following  sections  were  obtained.  The  grid 
systems  are  fully  3-D  (planar  symmetry  is  not  assumed),  model  both  the 
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Striker  assembly  as  well  as  the  catch  grill,  and  provide  an  acceptable 
balance  between  solution  accuracy  and  total  simulation  run  times.  As  will 
be  shown  in  more  detail  in  the  following  section,  the  final  grid  systems 
consisted  of  all  tetrahedral  cells  with  attention  given  to  gridding  around 
the  BLS  vent  systems,  the  striker  and  grill,  and  finally,  the  box  in  the  third 
case  involving  the  8-ftx8-ft  C2SQ  configuration.  A  large  plenum  was 
assumed  around  the  BLS  to  provide  a  space  into  which  the  air  flowing 
through  the  BLS  vents  can  exhaust  and  to  minimize  the  effects  of  waves 
reflecting  off  of  computational  outer  boundaries. 

Grids  were  generated  for  three  configurations  referred  to  in  this  appendix 
as  the  GSA  configuration,  8-ftx8-ft  C2SQ  configuration,  and  8-ftx8-ft 
C2SQ  configuration  with  box.  The  GSA  configuration  contained  all  stages 
of  the  BLS  configuration  through  the  circular  section.  The  second  and 
third  cases  added  the  square  8-ftx8-ft  test  section  to  the  GSA  configura¬ 
tion.  The  grid  for  the  third  case  was  modified  to  include  the  box  in  the 
C2SQ  section. 

All  analyses  were  run  on  the  MSU  High  Performance  Computing  Collabo- 
ratory  (HPCC)  IBM  IDataplex  and  Cray  CS300-LC  clusters  (named  Talon 
and  Shadow,  respectively).  Both  systems  use  multi-core  Intel  XEON  pro¬ 
cessors  for  compute  nodes.  A  complete  description  of  both  systems  is 
given  at  http:  / /hpc.msstate.edu.  Computational  grids  were  generated 
using  the  SolidMesh  grid  generation  tools  developed  by  Dr.  David  Marcum 
at  Mississippi  State  University  (Marcum  1998).  SolidMesh  provides  a  GUI- 
based  environment  that  supports  import  and  manipulation  of  CAD 
geometry  and  grid  generation  setup.  The  final  volume  grids  were 
generated  using  Marcum’s  Advancing  Front  with  Local  Reconnection  3-D 
(AFLR3-D)  tetrahedral  grid  generator.  All  calculations  discussed  herein 
were  made  assuming  inviscid  flow  and  used  the  ideal  gas  EOS.  Although 
Loci/CHEM  has  a  full  viscous  flow  analysis  capability,  finer  meshes  near 
solid  walls  would  be  required  to  resolve  boundary  layer  and  other  viscous 
effects,  which  in  turn,  would  require  switching  to  an  implicit  solver  instead 
of  the  Runge-Kutta  scheme  to  overcome  the  explicit  time-step  restrictions 
inherent  in  viscous  grids.  At  some  point  in  the  next  phase  of  this  research, 
an  attempt  should  be  made  to  quantify  viscous  effects,  but  time 
constraints  prohibited  any  viscous  analysis  in  the  current  simulations. 
Pressure-time-histories  were  extracted  at  the  physical  locations  of  the 
gauge  locations  prescribed  in  the  experiments.  These  were  integrated 
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using  trapezoidal  integration  in  a  MATLAB  program  to  obtain  time- 
histories  of  impulse  at  the  various  gauge  locations. 

Problem  setup 

Overview 

The  SolidMesh  program  was  used  to  define  the  surface  geometries  and 
meshes  of  the  three  BLS  configurations  along  with  computational  bound¬ 
aries  and  boundary  conditions.  Only  tetrahedral  meshes  were  used  to 
generate  the  3-D  results  presented  in  this  appendix.  The  final  tetrahedral 
volume  meshes  were  generated  using  the  AFLR3-D  program.  The  surface 
geometries  were  defined  using  a  series  of  drawings  and  dimensions  of  the 
various  BLS  sections  supplied  by  ERDC.  The  2-D  pathfinder  calculations 
used  both  triangular  and  quadrilateral  elements.  Initial  3-D  simulations 
for  the  first  case  did  not  include  either  the  striker  or  the  grill;  however, 
following  these  initial  simulations,  it  was  decided  to  include  these  in  all 
subsequent  3-D  calculations.  The  procedure  used  to  generate  the  meshes 
for  both  the  baseline  GSA  and  the  8-ftx8-ft  configurations,  as  well  as  the 
modeling  of  the  striker  mechanism  and  grill  are  described  in  the  following 
sections.  A  description  of  the  2-D  axisymmetric  and  initial  3-D 
calculations  for  both  the  GSA  and  8-ftx8-ft  configurations  are  given  in 
Mord  (2015).  For  this  report,  only  the  final  3-D  BLS  geometries  and 
meshes  are  discussed. 

Creation  of  baseline  meshes 

The  3-D  GSA  configuration  was  created  by  revolving  lines  about  the  x-axis 
to  form  a  set  of  surfaces  with  radial  symmetry.  The  8-ftx8-ft  configuration 
was  derived  from  the  GSA  model.  The  8-ftx8-ft  configuration  differs  from 
the  GSA  version  only  after  the  Cascade.  Therefore,  the  geometry  for  the 
symmetric  portion  of  the  BLS  inner  surface  was  created  from  a  set  of  lines 
revolved  around  the  x-axis.  The  non-symmetric  components  (C2sQ  and 
SQi)  were  developed  as  separate  surface  entities  and  merged  with  the 
revolved  line  components  to  form  the  final  BLS  model.  The  two  sections  of 
the  model  are  separated  by  a  small  gap  of  0.4  in.  This  distance  was  taken 
from  photographs  of  measurements  from  the  BLS  and  is  therefore 
consistent  with  the  BLS  experiments.  The  outlines  of  the  final  GSA  and  8- 
ftx8-ft  configuration  wall  geometries  along  the  plane  of  symmetry  can  be 
seen  in  Figures  Cl  and  C2. 
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Figure  Cl.  GSA  configuration  geometry. 


Figure  C2.  BLS  8-ftx8-ft  C2SQ  configuration  geometry  with  box. 


The  final  BLS  mesh  in  both  configurations  had  a  maximum  mesh 
resolution  of  i.o  in.  throughout.  The  extent  of  the  computational  grid  for 
the  8-ftx8-ft  configuration  with  box  is  seen  in  Figure  C3.  The  outer  mesh 
boundaries  were  extended  beyond  the  surface  of  the  BLS  to  provide  a  large 
region  into  which  the  air  moving  through  the  vents  can  exhaust  and  to 
place  the  computational  boundaries  at  a  sufficient  distance  from  the  BLS 
to  eliminate  any  potential  spurious  reflections  from  the  boundary  surfaces. 
The  point  spacing  was  reduced  to  0.5  in.  around  the  vents  and  Cascade  to 
accurately  capture  the  flow  around  the  edges  associated  with  each.  This 
spacing  was  chosen  because  the  walls  of  the  BLS  were  modeled  as  0.5-in. 
thick;  therefore,  this  was  the  coarsest  spacing  that  would  accurately 
capture  the  flow  around  the  edge.  For  the  final  meshes,  no  planar  or  axial 
symmetry  is  assumed  so  all  three  meshes  are  fully  3-D. 

Modeling  the  striker  mechanism  and  griii 

The  grill  and  striker  were  both  modeled  as  separate  entities  to  simplify 
grid  generation.  Modeling  these  as  separate  bodies  made  the  overall  grid 
generation  easier  with  the  added  bonus  that  both  the  grill  and  striker 
could  be  added  or  removed  without  compromising  the  rest  of  the  BLS. 
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Figure  C3.  BLS  8-ftx8-ft  C2SQ  configuration  with  box, 
mesh  in  Z=0  piane. 


The  striker  geometry  was  modeled  using  a  variety  of  methods.  The  section 
with  arms  was  made  using  trimmed  circular  surfaces.  The  needle  was 
made  as  a  thin  cylinder  without  a  point.  The  point  could  have  caused 
trouble  with  the  grid  size,  so  it  was  not  modeled.  The  spring  mechanism  on 
the  other  side  of  the  striker  was  also  modeled  as  a  solid  cylinder  due  to  the 
complexity  of  making  a  helical  spring.  The  point  spacing  around  the 
striker  was  defined  as  0.125  in.  because  a  small  distance  between  the 
striker  arms  and  the  wall  of  the  CYC  section  was  required,  and  0.125  in. 
was  already  the  required  spacing  of  the  grill;  theoretically,  a  point  spacing 
of  0.125  in.  was  the  smallest  value  that  could  be  used  without  significantly 
slowing  down  the  calculation.  The  mesh  around  the  striker  along  with  the 
driver,  grill,  and  initial  vents  along  the  plane  of  symmetry  can  be  seen  in 
Figure  C4.  A  close-up  view  of  the  mesh  around  the  striker  is  shown  in 
Figure  C5.  Figure  C6  shows  the  mesh  at  the  striker  face  looking  upstream. 

The  grill  surface  geometry  was  developed  by  first  creating  the  front  face  of 
the  grill  as  a  composite  surface  composed  of  many  smaller  surfaces.  This 
group  of  small  surfaces  was  then  copied  to  form  the  back  face  of  the  grill. 
Finally,  all  the  surfaces  that  connect  the  two  faces  were  defined.  The  spac¬ 
ing  around  the  grill  was  restricted  to  the  0.125-in.  thickness  of  the  grill 
walls.  This  is  similar  to  how  the  spacing  around  the  vents  was  defined;  the 
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flow  will  move  past  the  forward  and  backward  faces  of  the  grill  in  one  cell. 
Also,  the  grill  was  spaced  one  cell  width  away  from  the  BLS  inner  wall 
because  of  issues  that  occurred  with  attempts  to  physically  attach  the  grill 
to  the  wall.  The  mesh  around  the  grill  can  be  seen  in  Figure  C7. 


Figure  C4.  Driver  section,  computationai  grid  showing 
striker,  griii,  and  vents. 


Figure  C5.  Ciose-up  of  computationai  grid  around 
the  striker  mechanism. 


ERDC/GSL  TR-16-31 


118 


Figure  C6.  Cross  section  of  computationai  grid 
around  the  striker  mechanism. 


Figure  C7.  Computationai  grid  around  griii 
iooking  toward  driver  section. 


Case  specific  mesh  modifications  and  sizes 

The  mesh  used  in  the  initial  GSA  configuration  calculations  was  reduced 
from  101  million  to  52.4  million  cells  for  the  final  run.  It  included  the  addi¬ 
tion  of  the  models  for  the  striker  and  grill  that  were  neglected  in  the  origi¬ 
nal  mesh.  The  side-on  pressure  gauge  and  its  mounting  rod  were  not 
modeled.  Instead,  a  Loci/BLAST  data  probe  location  was  specified  that 
corresponded  to  the  location  of  the  side-on  gauge. 
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The  computational  mesh  for  the  8-ftx8-ft  C2SQ  configuration  (no  box) 
consisted  of  73.1  million  cells.  The  mesh  was  generated  using  the  geometry 
initially  developed  for  the  8-ftx8-ft  configuration  with  box  and  had  further 
refinements  around  the  striker  and  catch  grill.  An  attempt  was  made  to 
model  the  side-on  gauge  assembly  used  in  the  8-ftx8-ft  C2SQ  configu¬ 
ration  without  box  tests;  however,  the  mounting  rod  for  the  gauge  was  not 
included  due  to  problems  with  modeling  the  attachment  of  the  rod  to  the 
disk  and  the  entire  assembly  to  the  BLS  wall.  Instead,  the  side-on  gauge 
disk  was  fixed  in  space  and  floats  in  the  flow  path.  This  helped  reduce  the 
complexity  of  the  mesh  in  the  region  around  the  gauge  and  also  reduced 
the  overall  size  of  the  mesh  at  the  expense  of  some  uncertainty  in  the 
results  due  to  neglecting  the  mounting  rod. 

For  the  third  case,  the  box  was  formed  by  first  creating  the  five  wetted 
surfaces  of  the  box  and  then  trimming  the  four  vertical  faces  of  the  box 
with  the  BLS  floor.  The  point  spacing  around  the  target  box  was  chosen  to 
be  0.5  in.  A  close-up  view  of  the  mesh  around  the  box  is  shown  in 
Figure  C8.  This  mesh  contained  approximately  69.9  million  cells. 

Initial  conditions 

An  initial  absolute  pressure  of  1312.7  psi  (1298  psig)  was  specified  for  the 
driver  section.  Driver  and  ambient  temperatures  were  set  to  77°F.  Ambient 
pressure  in  all  other  regions  of  the  mesh  was  set  to  14.7  psi.  Initial  values 
of  density  were  computed  from  the  specified  pressure  and  temperature 
using  an  ideal  gas  EOS. 

Boundary  conditions 

A  second-order  reflecting  boundary  condition  in  which  values  in  ghost 
cells  used  to  implement  the  boundary  condition  are  obtained  by  reflecting 
adjacent  interior  values  is  imposed  on  the  exterior  boundary  surfaces,  the 
BLS  walls,  the  striker  assembly,  and  the  catch  grill. 
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Figure  C8.  Computational  grid  along  plane 
of  symmetry,  around  box. 


GSA  configuration,  caiibrations  Tests  30a  and  30b 

The  GSA  configuration  simulations  were  run  on  the  MSU  HPCC  Shadow 
cluster  using  200  processors.  The  full  simulation  consisted  of  several  runs 
with  restarts  and  ran  to  a  simulation  time  of  about  70  msec.  Total  wall 
clock  time  for  the  simulation  was  156.53  hr.  Loci/BLAST  computes  a  gross 
estimate  of  parallel  efficiency,  which  averaged  around  88  percent  for  the 
GSA  configuration  simulations  on  Shadow. 

The  results  for  this  case  were  first  analyzed  to  verify  that  a  planar  shock 
front  was  formed  in  the  test  section.  Pressure  contours  along  the  plane  of 
symmetry  and  the  magnitude  of  the  gradient  of  the  density  field  that 
provides  Schlieren-like  images  were  used  to  visualize  the  formation,  the 
movement  of  the  planar  shock  wave,  and  the  complex  flow  field  through 
the  driver,  striker,  grill,  and  CVC  sections.  These  images  were  generated 
using  the  Ensight  3-D  visualization  tool.  Pressure  contours  normalized  by 
the  maximum  reflected  pressure  at  the  end  plate  at  select  times  are  shown 
in  Figure  C9.  At  5  msec,  the  shock  front  has  exited  the  CVC  section  and  is 
about  one-quarter  the  distance  into  the  conical  section. 
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Figure  C9.  GSA  configuration  pressure  contours. 


(a)  5  msec 


(d)  25  msec 


The  shape  of  the  shock  front  is  slightly  bulged  along  the  centerline  and  is 
normal  to  the  BLS  wall.  High-pressure  gas  is  still  moving  through  the 
striker  and  grill.  The  flow  field  through  the  vent  system  is  extremely  com¬ 
plex  with  multiple  oblique  shocks  being  formed.  At  15  msec,  the  shock 
front  is  about  three-quarters  the  way  down  the  constant  area  section.  A 
planar  shock  front  is  clearly  seen  at  this  time.  At  about  20  msec,  the  shock 
front  has  reached  the  target  wall  and  is  reflected  back.  At  25  msec,  the 
reflected  shock  is  moving  upstream  and  has  reentered  the  constant  area 
section.  The  reflected  shock  continues  to  move  upstream  until  it  interacts 
with  the  residual  driver  gases. 

The  Schlieren-like  images  shown  in  Figure  Cio  give  an  alternative  view  of 
the  complexity  of  the  flow  field  during  all  phases  of  the  calculation.  In  par¬ 
ticular,  the  interaction  of  the  oblique  shocks  generated  by  the  flow  through 
the  vents  is  clearly  seen.  At  5  msec,  the  oblique  shocks  around  the  vents 
are  clearly  distinguished.  A  planar  shock  in  the  constant  area  section  is 
seen  at  15  msec.  The  oblique  shocks  from  the  vents  are  merging  into 
barrel-like  shock  structures  that  are  moving  into  the  conical  section.  At 
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Figure  CIO.  Density  gradient  magnitude 
contours,  GSA  configuration. 


(a)  5  msec 


(b)  15  msec 


(d)  30  msec 


25  msec,  the  reflected  shock  is  seen  to  be  moving  upstream  and  is  about 
to  interact  with  the  high-pressure  driver  gas  moving  downstream.  At 
30  msec,  the  reflected  shock  has  interacted  with  the  driver  gas  and  lost  all 
structure.  The  resulting  flow  field  at  this  time  is  highly  complex  and  most 
likely  consists  of  areas  of  highly  turbulent  flow  with  localized  vortices. 
Since  this  simulation  is  inviscid,  the  turbulent  nature  of  the  flow  field  and 
any  shock-boundary  layer  interactions  are  not  captured.  Therefore,  the 
actual  flow  field  will  be  different. 

Next,  the  Loci/BLAST  results  were  compared  with  the  measured  data  at 
the  P3,  P4,  and  PE  gauge  locations  (see  Figure  Cii).  The  computed  results 
were  time  shifted  by  -3.4  msec  to  match  the  measured  arrival  time.  As 
seen  in  the  figures,  Loci/BLAST  tends  to  overpredict  the  initial  peak  pres¬ 
sure.  There  are  also  small  differences  after  50  msec.  These  differences 
show  up  in  the  computed  impulse,  where  the  second  peak  is  slightly  over¬ 
predicted.  Even  with  these  differences,  one  notes  good  overall  agreement 
with  the  measured  data.  The  differences  are  probably  due  to  variations  in 
mesh  densitites  in  the  areas  around  the  gauge  locations. 
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Figure  Cll.  Time-shifted  gauge  comparisons  for 
Loci/BLAST,  GSA  configuration. 


(a)  Gauge  P3  (b)  Gauge  P4 


(c)  Gauge  PE 


Calibration  tests  with  8-ftx8-ft  C2SQ  configuration 

The  Loci/BLAST  calculation  for  this  case  was  run  on  the  MSU  HPCC 
Shadow  Cluster  using  200  processors.  As  previously  stated,  two  additional 
geometric  components  were  added  to  the  GSA  configuration  for  this  case. 
The  first  is  the  addition  of  the  C2sQ  and  SQi  sections.  The  second  is  an 
approximate  model  of  the  side-on  pressure  gauge  assembly.  The  analysis 
of  this  case  consisted  of  a  series  of  restarts  that  were  terminated  at  a  simu¬ 
lation  time  of  about  70  msec.  Total  wall  clock  time  for  the  simulation  was 
194.94  hr.  The  average  parallel  efficiency  estimate  was  88  percent.  As  in 
the  GSA  simulation,  visualization  of  pressure  contours  and  density  gradi¬ 
ent  magnitude  along  the  plane  of  symmetry  were  used  to  track  the  move¬ 
ment  of  the  planar  shock  through  the  BLS  and  to  examine  shock 
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interaction  with  the  driver  gas.  Pressure  contours  at  select  times  are 
shown  in  Figure  C12.  Results  at  previous  times  were  identical  to  the  GSA 
simulation  results.  At  20  msec,  the  planar  shock  is  seen  to  have  moved 
through  the  Cascade  and  into  the  C2SQ  transition  section.  The  shock  wave 
is  approaching  the  back  end  of  the  SQi  section  at  25  msec.  The  reflected 
shock  wave  is  moving  upstream  at  35  msec  and  has  exited  the  C2SQ 
transition  section.  At  40  msec,  the  reflected  wave  has  moved  into  the  GSA 
section  and  is  interacting  with  the  driver  gas  moving  downstream. 

Figure  C12. 8  ftx8  ft  C2SQ  configuration 
pressure  contours. 


31 


(a)  20  msec 


(b)  25  msec 


(c)  35  msec 


(d)  40  msec 

Schlieren-like  images  are  shown  at  select  times  in  Figure  C13.  As  with  the 
first  case,  these  images  illustrate  the  complexity  of  the  flow  field  particu¬ 
larly  when  the  reflected  shock  interacts  with  the  driver  gas.  The  incident 
shock  wave  at  20  msec  and  the  reflected  wave  at  35  msec  are  clearly  seen. 
At  40  msec,  the  reflected  wave  is  being  eaten  away  by  the  driver  gases 
along  the  BLS  centerline.  The  existence  of  density  gradients  downstream 
of  the  interaction  again  indicates  the  complexity  of  the  flow  field  in  this 
region. 
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Figure  C13.  Density  gradient  magnitude  contours, 
8-ftx8-ft  C2SQ  configuration. 


(c)  40  msec 


The  calculated  pressure-time-histories  and  corresponding  impulse  are 
compared  with  measured  data  in  Figure  C14.  The  computed  results  were 
time  shifted  by  -3.6  msec  to  match  the  measured  arrival  time.  Compari¬ 
sons  are  made  for  the  Pi,  P5,  P6,  and  PR  (side-on)  pressure  gauge  loca¬ 
tions.  The  Loci/BLAST  results  for  gauges  Pi,  P5,  and  P6  are  in  excellent 
agreement  with  the  measured  data.  The  initial  peak  pressure  and  rise  time 
are  very  close  to  that  measured.  The  computed  waveforms  at  these  gauges 
are  almost  identical  to  the  experimental  results  except  for  small  variations 
in  the  pressure  oscillations.  The  integrated  impulse  also  compares  well. 
Similar  results  are  also  obtained  for  the  side-on  gauge  PR.  The  biggest 
anomoly  is  the  large  overprediction  of  the  pressure  spike  around  40  msec. 
Other  than  this  anomoly,  Loci/BLAST  is  seen  to  provide  an  accurate  pre¬ 
diction  of  the  waveform  for  the  side-on  gauge. 
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Figure  C14.  Time-shifted  gauge  comparisons,  8-ftx8-ft 
C2SQ  configuration. 


(a)  Gauge  PI  (b)  Gauge  P5 


(c)  Gauge  P6  (d)  Gauge  PR 

BLS  8-ftx8-ft  C2SQ  configuration  with  box 

The  Loci/BLAST  calculation  for  this  case  was  started  on  the  MSU  HPCC 
Talon  system.  When  the  Shadow  system  became  available,  it  was  restarted 
on  Shadow  and  run  out  to  a  simulation  time  of  around  65  msec.  The  cal¬ 
culation  ran  for  47.5  hr  on  Talon  using  192  processors  and  214.2  hr  on 
Shadow  using  200  processors  for  a  total  run  time  of  261.7  hr.  However,  the 
XEON  processors  used  on  Shadow  are  at  least  twice  as  fast  (probably 
more)  as  the  processors  used  on  Talon,  and  the  time  required  to  run  the 
entire  problem  on  Shadow  is  estimated  to  be  about  238  to  240  hr  of  wall 
clock  time.  The  average  parallel  efficiency  was  about  88  percent  on 
Shadow. 


Impulse  (psi-msec)  Impulse  (psi-msec) 
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As  with  the  other  test  cases,  pressure  contours  and  density  gradient  mag¬ 
nitude  contours  were  plotted  to  track  the  shock  movement  through  the 
BLS  and  to  visualize  the  complex  flow  field  around  the  box  when  both  the 
incident  and  reflected  shock  waves  intersect  and  pass  over  it.  The  pressure 
contours  are  for  pressures  normalized  by  an  estimate  of  the  peak  reflected 
pressure  on  the  face  of  the  target  at  pressure  gauge  PF3.  Figure  C15  shows 
the  pressure  contours  at  select  times.  Comparison  with  the  results  from 
the  second  case  suggests  that  the  box  had  little  effect  on  the  location  of  the 
incident  and  reflected  shock  at  the  same  simulation  times. 


Figure  C15.  Pressure  contours,  8-ftx8-ft 
C2SQ  configuration  with  box. 


(a)  20  msec 


(c)  35  msec 


(d)  40  msec 


The  Schlieren-like  images  for  this  case  showed  similar  results  for  the  inci¬ 
dent  and  reflected  waves  as  seen  in  the  second  case.  In  order  to  get  a  more 
detailed  view  of  the  interaction  of  the  incident  and  reflected  shock  waves 
with  the  box,  close-up  views  of  the  pressure  contours  were  made  for  times 
during  which  the  incident  and  reflected  shocks  passed  over  the  box.  The 
images  for  the  incident  shock  wave  are  shown  in  Figure  C16.  The  pressure 
in  the  figure  has  been  normalized,  based  on  the  peak  pressure  on  the  front 
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Figure  C16.  Contours  of  incident  shock  pressure  around  box. 


(a)  21  msec 


(b)  21.5  msec 


(c)  22  msec  (d)  22.5  msec 

face  (specifically  with  respect  to  gauge  PBF3).  The  passage  of  the  shock 
over  the  box  along  with  the  reflection  off  of  the  front  face  and  reformation 
of  the  shock  after  it  passes  the  back  face  of  the  box  are  clearly  seen.  The 
shock  front  moves  around  the  box  reaching  the  top  edge  of  the  back  face 
around  22  msec.  After  passage,  the  bottom  half  of  the  shock  lags  the  upper 
part  with  a  noticable  bend  around  the  upper  back  edge  of  the  box.  The 
reflected  wave  from  the  front  face  is  seen  to  rapidly  move  upstream  and 
begin  to  lose  strength  as  it  moves  toward  the  driver. 


Similar  images  were  generated  for  the  interaction  of  the  primary  reflected 
wave  with  the  box.  Figure  C17  shows  the  pressure  contours  in  the  vicinity 
of  the  box  at  select  times.  The  pressure  in  the  figure  has  been  normalized 
based  on  the  peak  pressure  on  the  front  face  (specifically  with  respect  to 
gauge  PBF3).  As  with  the  incident  shock  interaction,  the  primary  reflected 
shock  also  produces  a  secondary  reflected  shock  that  moves  downstream 
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(a)  31.5  msec 


(c)  32.5  msec 


(b)  32  msec 


(d)  33  msec 


(e)  33.5  msec 

from  the  back  face  of  the  box.  As  the  primary  reflected  wave  moves  over 
the  box  it  assumes  a  bow-like  shape  due  to  interaction  with  the  box.  It 
eventually  returns  to  a  more  planar  shape  as  it  moves  upstream  out  of  the 
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transition  region  until  it  begins  to  interact  with  the  driver  gas,  after  which 
it  loses  a  coherent  shape. 

Computed  pressure-time-histories  and  resulting  impulse  were  compared 
against  the  measured  data  at  several  gauge  locations  on  the  front,  back, 
top,  left,  and  right  side  faces  of  the  box.  In  addition,  comparisons  were 
drawn  against  P3  and  P9  gauges  located  on  the  target  wall.  The  computed 
results  were  time  shifted  by  -3.1  msec  to  match  the  measured  arrival  time. 
The  results  for  front  face  gauges  are  shown  in  Figure  C18.  The  correlation 
with  the  measured  data  are  excellent.  The  computed  waveform  shapes  are 
similar  to  that  measured  for  all  gauges  except  PBF6,  where  Loci/  BLAST  is 
seen  to  slightly  overpredict  the  incident  peak  pressure.  Overall,  the 
secondary  shock  at  30  msec  is  accurately  captured  with  only  minor 
differneces  noted  after  50  msec. 

The  results  for  the  back  face  gauges  are  shown  in  Figure  C19.  The  correla¬ 
tion  with  the  measured  data  are  also  excellent.  As  with  the  front  face 
gauges,  the  overall  waveforms  are  accurately  predicted  with  the  largest 
variation  seen  in  the  results  for  PBB2  and  PBB4  where  the  peak  pressure  is 
slightly  overpredicted;  however,  this  has  little  effect  on  the  integrated 
impulse. 

The  results  for  the  left  face  gauges  are  shown  in  Figure  C20.  The  excellent 
correlation  with  experimental  data  seen  in  the  front  and  back  face  results 
are  repeated  for  PBLi  and  PBL2.  There  are  no  comparisons  for  PBL3  and 
PBL4.  Following  the  analysis,  it  was  found  that  these  gauges  were 
misplaced  in  the  model.  The  setup  assumed  symmetry  of  the  gauge 
positioning  for  the  right  and  left  sides;  however,  this  was  not  the  case  in 
the  experimental  setup  where  the  third  and  fourth  gauges  were 
asymmetric  for  the  two  sides. 

The  results  for  the  right  face  gauges  are  shown  in  Figure  C21.  All  of  these 
gauges  were  correctly  placed  in  the  model  and  match  the  locations  in  the 
experimental  setup.  Once  again,  the  computed  waveforms  are  similar  to 
the  measured  ones. 

Figure  C22  provides  a  comparison  for  the  gauges  on  the  top  face.  Once 
again,  good  correlation  is  noted  with  the  measured  data.  Peak  pressures 
are  extremely  close  for  all  gauges  along  with  the  impulse. 
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Figure  CIS.  Time-shifted  gauge  comparisons 
for  Loci/BLAST,  front  face. 
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Figure  C19.  Time-shifted  gauge  comparisons 
for  Loci/BLAST,  back  face. 
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Figure  C20.  Time-shifted  gauge  comparisons 
for  Loci/BLAST,  ieft  face. 
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Figure  C21.  Time-shifted  gauge  comparisons 
for  Loci/BLAST,  right  face. 
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Figure  C22.Time-shifted  gauge  comparisons 
for  Loci/BLAST,  top  face. 
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Finally,  the  comparisons  for  gauges  P3  and  P9  are  shown  in  Figure  C23. 
The  computed  results  for  the  P3  gauge  overpredict  the  peak  pressure  but 
still  captures  the  measured  waveform  with  good  accuracy.  The  results  for 
the  P9  gauge  are  in  much  better  agreement  than  the  P3  gauge.  Again,  the 
measured  waveform  shape  is  accurately  reproduced. 

Figure  C23.  Time-shifted  gauge  comparisons 
for  Loci/BLAST,  target  waii. 
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Conclusions  and  recommendations 

The  MSU  Loci/BLAST  code  was  used  to  successfully  simulate  the  3-D  flow 
fields  in  the  ERDC  BLS  for  three  different  experimental  configurations. 
Correlations  were  made  with  measured  pressure-time-histories  and 
impulse  values  at  a  variety  of  gauge  locations.  In  general,  the  quality  of  the 
comparisons  ranged  from  very  good  to  excellent.  In  particular,  the  results 
for  the  8-ftx8-ft  C2SQ  configuration  with  box  demonstrated  that  Loci/ 
BLAST  can  accurately  reproduce  both  incident  and  reflected  planar  wave 
interactions  with  the  box.  In  most  cases,  the  calculated  waveforms  were 
similar  to  the  measured  ones.  These  results  indicate  that  Loci/BLAST  is  an 
excellent  candidate  for  future  BLS  analyses  and  related  airblast  analyses. 

Based  on  these  results,  the  following  recommendations  have  been  made. 
The  grid  size  used  in  the  current  simulations  varied  from  about  50  to 
70  million  cells.  A  much  larger  mesh  needs  to  be  generated  and  run  on  the 
ERDC  HPC  systems  to  quantify  the  improvements  in  accuracy  that  can  be 
obtained  by  a  much  finer  mesh.  This  will  also  allow  us  to  investigate  the 
parallel  scalability  of  Loci/BLAST  and  Loci/CHEM  for  both  BLS  analyses 
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and  large-scale  airblast  problems.  Since  all  the  simulations  in  the  current 
research  were  inviscid  and  neglected  important  viscous  effects  such  as 
turbulence  and  shock-boundary  layer  interactions,  a  viscous  analysis 
should  be  attempted  in  the  next  phase  of  this  project.  The  excellent  corre¬ 
lations  obtained  with  the  inviscid  solutions  indicate  that  these  effects 
might  be  secondary  for  the  current  test  cases.  However,  viscous  effects  still 
need  to  be  examined  to  quantify  under  what  circumstances  they  become 
important. 
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Appendix  D;  RAGE  Modeling 

RAGE  background 

The  Radiation  Adaptive  Grid  Eulerian  (RAGE)  code  was  developed  by 
Leidos,  formerly  Science  Applications  International  Corporation  (SAIC) 
(Gittings  et  al.  2008).  There  was  no  user  manual  available  for  RAGE;  how¬ 
ever,  there  was  a  draft  manual  and  example  problem  set  available  for  a 
companion  code,  the  SAIC  Adaptive  Grid  Eulerian  (SAGE)  (Archer  et  al. 
2005).  SAGE  has  been  in  use  at  ERDC  for  many  years,  and  there  is  consid¬ 
erable  in-house  experience  modeling  high  explosive  events  with  this  code. 
This  provided  a  launching  point  for  running  RAGE.  Fortunately,  both 
SAGE  and  RAGE  shared  a  common  set  of  input  commands  and  structure. 
This  made  the  transition  from  running  SAGE  to  RAGE  fairly  seamless.  The 
ensuing  discussion  of  RAGE  draws  upon  that  in  Gittings  et  al.  (2008)  and 
Archer  et  al.  (2005),  as  well  as  our  past  experience  with  SAGE  at  ERDC. 

One  should  note  that  there  are  actually  two  versions  of  RAGE,  one  main¬ 
tained  by  Leidos  and  the  other  by  Los  Alamos  National  Laboratory  (LANL) 
(Clover  2015).  The  Leidos  version  of  RAGE  was  utilized  for  all  calculations 
discussed  herein.  The  development  path  for  the  code  base  appears  to  have 
split  in  the  2006  time  frame.  All  calculations  discussed  herein  were 
conducted  on  the  ERDC  HPC  platform  Garnet  using  Version 
20110102.000,  which  was  installed  by  Leidos  at  the  onset  of  this  project. 
Garnet  is  a  Cray  XE6  that  has  4,716  compute  nodes  with  32  cores  per 
node.  Each  compute  node  has  a  core  speed  of  2.5  GHz  and  64  GB  of 
accessible  memory. 

RAGE  is  a  multidimensional,  massively  parallel  code  that  uses  a  two-shock 
or  single-intermediate-state  approximate  Riemann  solver  in  the  shock 
capturing.  Both  minmod  and  modified  van  Leer  slope  limiters  are 
available  for  use  in  calculating  gradients  across  cell  faces.  The  latter  sets 
the  flux  limiter  to  1.5  and  has  the  benefit  of  reducing  ringing  behind  the 
shock  front.  RAGE  has  the  ability  to  handle  multi-material  cells  and  a 
robust  adaptive  mesh  refinement  (AMR)  capability.  The  AMR  capability 
was  of  great  interest  in  this  work,  as  a  fine  resolution  is  desired  when 
modeling  shock  propagation. 

All  RAGE  calculations  are  conducted  on  a  structured  block  mesh  com¬ 
posed  of  rectangular  cells,  i.e.,  a  mesh  of  fixed  spacing  along  the  Cartesian 
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axes.  The  number  of  cells  assigned  to  a  particular  axis  can  vary;  however, 
the  cell  count  must  be  even  to  support  AMR.  The  code  supports  i-D,  2-D, 
and  3-D  problem  setups,  all  within  the  same  architecture.  It  uses  an 
h- adaptive  approach  in  the  mesh  refinement,  where  cells  are  split  and/or 
joined  uniformly  based  on  an  indicator.  All  mesh  adaption  is  conducted 
subject  to  a  2-to-i  refinement  rule.  This  rule  requires  that  the  level  of 
refinement  not  differ  by  more  than  one  level  between  adjacent  cells.  This 
rule  applies  to  both  refinement  and  unrefinement.  The  mesh  adaption  is 
continuous  in  both  space  and  time,  with  cell  adaption  triggered  based  on 
changes  in  state  data,  e.g.,  velocity  or  pressure,  between  cells  or  via  esti¬ 
mates  of  the  first-order  truncation  error  (Gittings  et  al.  2008).  Details  on 
the  indicator  thresholds  were  not  available.  Based  on  our  personal  experi¬ 
ence  with  the  code,  it  appears  the  refinement  strategy  is  geared  towards 
shock-following. 

For  structured  Eulerian  codes,  boundary  conditions  are  applied  to  the 
mesh  extents,  rather  than  individual  materials.  The  available  documenta¬ 
tion  indicated  that  only  reflective  boundary  conditions  are  supported. 
Further  discussion  with  Leidos  indicated  that  outflow  boundary  conditions 
were  available  as  well;  however,  these  were  not  suited  to  the  particular 
setup  chosen  for  the  BLS,  which  involved  the  use  of  rigid  materials  (Clover 
2015).  Thus,  the  analysis  setup  relied  on  reflective  boundaries  at  the  mesh 
extents.  This  boundary  condition  type  can  lead  to  artificial  reflections  from 
the  mesh  extents  that  can  affect  the  solution  later  in  time.  The  only  way  to 
avoid  this  is  to  extend  the  mesh  well  beyond  the  physical  bounds  of  the 
BLS.  This  helps  mitigate  any  artificial  reflections  at  the  cost  of  additional 
computations  on  the  extended  mesh. 

Materials  are  inserted  into  the  mesh  to  develop  the  problem  setup.  This 
generally  takes  the  form  of  overlaying  primitive  objects  (blocks,  cones,  and 
spheres  for  3-D  problems)  to  build  up  a  desired  geometry.  The  order  of  the 
overlay  process  determines  whether  a  material  is  inserted  (added)  or 
masked  (cutout).  The  walls  of  the  BLS  were  modeled  as  rigid  (non¬ 
responding)  materials.  The  surrounding  air  was  modeled  as  a  hydrody¬ 
namic  material.  It  should  be  noted  that  RAGE  supports  geometry  import 
from  a  CAD  package;  however,  this  capability  was  not  utilized  for  the 
calculations  discussed  herein. 

The  Leidos  version  of  RAGE  used  in  these  analyses  does  not  have  a  mate¬ 
rial  interface  reconstruction  algorithm.  The  consequence  is  that  increased 
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resolution  may  be  needed  to  ensure  that  the  thin  shell  of  the  BLS  structure 
is  captured  in  the  modeling.  Meeting  this  resolution  requirement  can  lead 
to  excessively  long  run  times  that  are  intractable  for  3-D  analyses.  It  is  not 
just  a  matter  of  resolving  the  shell  thickness,  but  also  the  need  to  resolve 
the  shell  along  the  full  length  of  the  BLS  (about  40  ft)  as  well  as  shocks 
traversing  the  extended  mesh  outside  the  BLS.  The  AMR  capability  is 
helpful  to  efficiently  partition  the  refined  regions  initially;  however,  as  the 
calculation  progresses,  the  number  and  complexity  of  the  shock  reflections 
increase  resulting  in  growth  of  the  overall  problem  size.  It  is  up  to  the 
analyst  to  achieve  a  balance  between  the  desired  mesh  resolution, 
available  computational  resources,  and  turnaround  time  for  the  analysis. 

There  are  two  approaches  for  dealing  with  the  issue  of  capturing  the  thin 
shell  of  the  BLS  with  an  acceptable  resolution.  The  first  is  to  modify  the 
location  of  the  shell  structure  so  that  the  material  interfaces  are  aligned 
with  the  cell  boundaries.  The  second  is  to  artificially  increase  the  thickness 
of  the  shell  structure  to  allow  more  cells  to  span  its  thickness  while  at  a 
lower  level  of  resolution.  Both  options  were  investigated;  however,  one 
should  note  that  the  ability  to  carry  out  the  former  is  difficult  given  the 
conical  nature  of  several  BLS  sections. 

RAGE  has  several  EOS  modeling  options  suitable  for  airblast  applications. 
These  include  SESAME  tables  for  air  and  an  extensive  EOS  library  that 
covers  a  range  of  conventional  explosives.  The  options  for  modeling  the 
explosive  process  are  limited  to  an  ignition  and  growth  model  (Lee  and 
Tarver  1980),  which  transitions  the  explosive  from  its  initially  unreacted 
to  fully  reacted  states  both  described  using  a  Jones- Wilkins-Lee  (JWL) 
EOS.  RAGE  also  supports  an  ideal  gas  EOS.  Both  the  SESAME  and  ideal 
gas  EOS  options  were  of  interest  for  the  BLS  modeling. 

RAGE  requires  that  all  data  be  in  the  centimeter-gram-second  (CGS)  sys¬ 
tem  of  units.  With  this  system  of  units,  pressure  is  in  units  of  dynes/cm^, 
density  in  g/cm3,  velocity  in  cm/s,  and  temperature  in  electron-volts.  For 
presentation  purposes,  the  RAGE  results  are  typically  cast  in  English  units 
to  maintain  consistency  with  the  overall  report  write-up.  There  are  select 
circumstances  where  the  data  are  presented  in  CGS  units,  mainly  due  to 
limitations  in  modifying  plot  ranges.  The  units  will  be  clearly  defined  when 
this  occurs. 
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GSA  configuration,  caiibration  Tests  30a  and  30b 

The  initial  RAGE  modeling  involved  a  series  of  2-D  cylindrical  (2-DC) 
calculations.  The  first  set  of  calculations  involved  a  resolution  study,  while 
the  second  looked  at  the  sensitivity  of  various  modeling  options  on  the 
predictions.  The  problem  setup  used  for  the  2-DC  calculations  is  shown  in 
Figure  Di.  The  image  represents  only  a  small  portion  of  the  computational 
grid.  Only  the  axisymmetric  components  of  the  BLS  can  truly  be  modeled 
in  a  2-DC  calculation.  Consequently,  neither  the  grill  nor  leg  supports  in 
the  striker  flange  can  be  properly  modeled.  It  was  felt  that  the  opening  (or 
vent)  area  in  the  striker  flange  should  be  preserved.  To  do  so,  the  diameter 
of  the  ring  support  at  the  base  of  the  striker  mechanism  was  artificially 
increased  to  maintain  an  equivalent  opening  area. 


Figure  Dl.  RAGE  model  setup  for  2-DC  calculations 
involving  the  GSA  configuration. 
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The  computational  grid  extended  well  beyond  the  BLS  to  minimize  any 
effects  from  artificial  reflections  at  the  mesh  boundaries.  The  initial  cell 
size  ho  was  set  to  64  cm  for  all  calculations  discussed  herein.  When  invok¬ 
ing  the  AMR  capability,  the  user  must  also  specify  the  maximum  level  of 
resolution  (MLOR).  Once  these  input  data  are  specified,  one  can  then 
determine  the  size  of  the  finest  cell  in  the  mesh  hnwR  as  follows. 

h-MLOR  ~  2(MLOR-i)  (D"1) 

For  the  2-DC  calculations,  the  mesh  extended  out  40  cells  in  the  axial  and 
lateral  directions  with  the  origin  defined  at  the  base  of  the  driver.  Thus, 
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the  mesh  encompassed  a  physical  space  that  was  2,560  cm  by  2,560  cm 
(40  times  ho  in  either  direction),  with  the  longitudinal  axis  of  the  BLS 
aligned  with  the  axis  of  revolution  in  the  problem  setup.  The  extended 
mesh  was  thought  sufficient  to  avoid  any  reflections  that  might  originate 
from  the  mesh  boundaries  and  affect  the  solution  later  in  time. 

The  walls  of  the  BLS  were  modeled  as  rigid  (non-responding)  materials. 
The  outer  shell  thicknesses  of  the  CVC  sections  depicted  in  Figure  Di  were 
artificially  increased  to  a  thickness  of  0.79  in.,  about  1.6  times  the  true 
thickness.  A  trial-and-error  approach  was  taken  in  determining  the 
increased  shell  thickness.  A  set  of  scoping  calculations  was  conducted 
using  a  MLOR  of  9  (0.25-cm  cell  size)  and  the  true  shell  thickness.  The 
thickness  was  increased  in  0.1-in.  increments,  until  it  was  determined  that 
the  shell  structure  was  adequately  captured  in  the  model.  The  choice  of  a 
starting  MLOR  of  9  was  subjective,  but  seemed  reasonable  based  on  expe¬ 
rience  with  other  Eulerian  codes  and  their  requirements  for  resolution 
through  the  thickness  of  a  material.  Increasing  the  CVC  shell  thickness  to 
0.79-in.  results  in  a  decrease  in  vent  area;  this  is  on  average  7.5  percent  for 
each  vent.  This  reduction  should  have  no  effect  on  the  peak  pressures  for 
the  target  wall  gauges;  however,  it  could  potentially  affect  the  late-time 
blow-down  of  the  gas,  and  hence  the  impulse. 

The  shell  thickness  for  the  driver  and  other  downstream  sections  was 
increased  slightly  more  to  about  1.2-in.  thick.  The  artificially  increased 
thickness  for  these  other  sections  should  have  no  effect  on  the  calculations 
since  they  have  no  vents  and  only  serve  to  contain  the  gas  in  the  BLS. 

The  material  surrounding  the  BLS  structure  was  composed  of  air.  Only  a 
single  air  material  was  specified  throughout;  however,  the  initial  condi¬ 
tions  varied  by  region  to  match  those  specified  in  the  test.  Here,  the  pres¬ 
sure  of  the  driver  gas  was  set  to  1,298  psig  and  was  assumed  to  be  at  an 
ambient  temperature  of  77°F.  The  remainder  of  the  air  both  inside  and 
outside  of  the  BLS  was  set  to  ambient  conditions  (14.7  psia,  77°F).  The  air 
was  modeled  as  an  ideal  gas  having  a  ratio  of  specific  heats  (gamma  y)  of 
1.4,  which  is  reasonable  for  the  range  of  conditions  considered  here.  The 
air  for  all  calculations  was  modeled  as  an  ideal  gas.  It  was  desired  to  also 
consider  a  SESAME  EOS  for  real  air;  however,  it  was  found  that  run  times 
with  this  EOS  were  intractable.  It  was  decided  that  all  modeling  would  rely 
solely  on  an  ideal  gas  EOS. 
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A  mesh  convergence  study  was  conducted.  The  results  are  shown  in 
Figure  D2  for  three  gauge  locations,  P3,  P4,  and  PE.  The  radial  distance 
from  the  BLS  centerline  for  these  gauges  was  2.25, 13.2,  and  30  in., 
respectively.  The  gauges  selected  provide  a  good  indication  of  the 
influence  of  radial  offset  on  the  data.  In  practice,  the  measured  (and 
calculated)  results  were  generally  consistent  across  the  target  wall.  The 
resolution  for  each  calculation  is  denoted  in  parentheses,  along  with  the 
MLOR.  Overall,  one  notes  that  the  calculated  TOA  is  slightly  larger  than 
measured  (on  the  order  of  20  percent)  with  the  initial  peak  pressure  being 
significantly  over-predicted  (on  the  order  of  35  to  40  percent).  In  turn,  this 
leads  to  a  notable  overprediction  in  the  impulse. 


Figure  D2.  Comparison  of  2-DC  RAGE  resuits  for  GSA  configuration. 
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The  calculated  results  at  gauge  P3  are  exceedingly  noisy.  This  gauge  is  near 
the  centerline  of  the  BLS  (2.25-in.  offset).  The  BLS  walls  and  any 
symmetry  planes  in  a  calculation  act  as  perfect  reflectors,  and  it  is  conjec¬ 
tured  that  this  idealization  induced  ringing  in  the  solution  due  to  the  axi- 
symmetric  nature  of  the  problem  setup.  The  ringing  in  the  calculated 
result  diminishes  with  increased  radial  distance  from  the  centerline.  Thus, 
the  results  for  gauges  P4  and  PE  provide  the  greatest  insight  on  the  influ¬ 
ence  of  mesh  resolution.  It  is  clear  that  only  marginal  improvements  in  the 
impulse  are  noted  with  increasing  resolution.  A  sub-millimeter  resolution 
is  certainly  more  than  adequate  to  capture  the  shock  propagation.  Thus, 
one  must  conclude  that  the  overprediction  is  due  to  some  other  factor  in 
the  problem  setup. 

A  sensitivity  study  was  conducted  next  to  further  investigate  the  problem 
setup.  The  results  of  the  study  are  outlined  for  gauge  PE  in  Figure  D3. 
These  results  are  typical  of  those  noted  at  the  other  gauge  locations.  All 
calculations  are  based  on  a  MLOR  of  9  (0.25-cm  cell  size).  They  utilized 
the  mock-up  striker  mechanism  depicted  in  Figure  Di,  unless  stated  other¬ 
wise.  Also,  the  modified  van  Leer  limiter  was  used  throughout,  unless 
stated  otherwise. 

The  first  variation  considered  the  influence  of  the  flux  limiter  (see 
Figure  D3(b)),  where  the  less  accurate  minmod  slope  limiter  was  specified. 
This  did  not  influence  the  calculated  results.  The  next  variation  did  not 
include  any  striker  mechanism  in  the  setup,  thus  focusing  on  the  influence 
of  the  striker  assembly  on  the  downstream  flow  field  (see  Figure  D3(c)). 
Here,  a  slight  increase  in  the  peak  reflected  pressure  was  noted  when  the 
striker  was  omitted  from  the  setup.  Recall,  the  grill  is  not  included  in  the 
2-DC  model.  Thus,  one  might  infer  that  the  grill  may  lead  to  a  constriction 
in  flow,  and  hence  a  reduction  in  the  gauge  pressures  measured  down¬ 
stream.  As  an  aside,  it  is  useful  to  review  an  earlier  test  that  did  not  have  a 
grill  (see  Figure  D4).  The  TOA  for  Test  24  in  the  figure  was  time-shifted  to 
match  that  for  Test  30a.  There  were  issues  triggering  the  measured  data  in 
the  earlier  testing.  Although  the  measured  TOA  from  these  earlier  tests  is 
suspect,  the  general  shape  of  the  pulse  is  not.  It  is  clear  that  the  two  sets  of 
data  compare  well,  suggesting  that  the  grill  does  not  significantly  influence 
the  measured  data.  One  can  then  infer  that  the  grill  would  not  affect  the 
calculations,  had  it  been  included  in  the  problem  setup. 


Overpressure  (psi)  Overpressure  (psi)  Overpressure  (psi) 
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Figure  D3.  Results  from  2-DC  sensitivity  study. 


(a)  All  variations 


(b)  Flux  limiter 


Time  (msec) 


Time  (msec) 


(c)  No  striker  mechanism 


Time  (msec) 


(d)  Separate  gas  and  ambient  air 


Time  (msec) 


(e)  Driver  pressure 


(f)  Ratio  of  specific  heats 
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Figure  D3.  (Continued). 


Time  (msec) 


(g)  Grid  aligned  BLS  structure 


Figure  D4.  Comparison  of  gauge 
PE  measured  data  with 
and  without  grill. 


Further  variations  considered  the  air  modeled  as  two  separate  materials, 
ambient  air  and  the  compressed  driver  gas  (see  Figure  D3(d));  variations 
in  the  initial  driver  pressure  as  ±2  atmospheres  (+29.4  psi)  as  shown  in 
Figure  D3(e);  variations  in  gamma  by  ±5  percent  as  shown  in  Figure  D3(f); 
and  modifying  the  location  of  the  BLS  outer  shell  structure  so  it  was 
aligned  on  cell  boundaries  (see  Figure  D3(g)).  Note  that  it  is  not  possible 
to  align  the  conical  segments  along  the  cell  boundaries  as  they  are  non- 
orthogonal  entities.  The  results  were  fairly  insensitive  for  all  cases  except 
for  variations  in  gamma. 
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One  should  note  that  the  ratio  of  specifics  heats,  or  gamma,  is  assumed 
constant  throughout  the  entire  domain,  regardless  of  the  changing 
conditions.  In  reality,  one  can  expect  some  variation  in  gamma;  however, 
it  should  be  reasonably  small  for  the  limited  range  of  conditions  consid¬ 
ered  here.  The  bounding  cases  for  gamma  used  in  the  sensitivity  study  are 
not  applicable  to  the  case  at  hand;  however,  they  do  indicate  that  even 
with  these  unrealistic  choices,  the  initial  peak  pressure  is  still  significantly 
overpredicted. 

There  are  other  possible  causes  for  the  differences  between  the  calculated 
and  measured  results.  These  include  modeling  (or  lack  of)  diaphragm 
rupture  and  the  reduced  venting  area  in  the  CVC  (necessary  to  conduct  the 
calculation).  The  driver  gas  is  assumed  to  be  at  a  uniform  pressure.  When 
the  calculation  starts,  the  driver  gas  instantaneously  flows  downstream 
with  no  impedance.  It  is  possible  that  the  process  of  rupture  and  lack  of 
uniformity  in  the  initial  flow  may  have  some  effect;  however,  it  is  doubtful 
that  the  effect  would  be  as  great  as  observed  in  these  calculations.  Thus, 
the  reasoning  for  the  overpredictions  is  at  present  unknown. 

Although  the  GSA  configuration  is  predominately  axisymmetric,  the 
emphasis  is  on  3-D  modeling  as  this  provides  a  better  indication  of  a 
code’s  performance  for  handling  realistic  applications.  Thus,  all  further 
calculations  involved  modeling  the  BLS  configurations  in  3-D.  Note, 
symmetry  conditions  about  the  longitudinal  axis  were  considered  to 
reduce  the  run  time.  The  2-DC  sensitivity  study  suggests  that  RAGE  will 
overpredict  both  the  peak  pressure  and  impulse  in  the  3-D  setup.  Given 
the  time  constraints  on  the  project,  it  was  necessary  to  proceed  forward 
with  the  3-D  calculations  with  the  expectation  that  both  the  peak  pressure 
and  impulse  would  be  overpredicted. 

The  3-D  problem  setup  for  the  GSA  configuration  is  outlined  in  Figure  D5. 
Quarter  symmetry  was  assumed  with  symmetry  planes  about  the  vertical 
and  horizontal  planes  along  the  longitudinal  BLS  axis.  The  figure  depicts 
only  a  portion  of  the  mesh.  The  actual  mesh  extended  out  along  the  cardi¬ 
nal  directions  (30  cells  along  the  x-and  y-axes  and  40  cells  along  the 
z-axis).  As  with  the  2-DC  setup,  the  origin  was  defined  at  the  base  of  the 
driver.  The  initial  cell  size  was  64  cm,  resulting  in  a  grid  encompassing  a 
1920-cm  X  1920- X  2560-cm  region  of  physical  space.  The  extended  mesh 
was  thought  sufficient  to  avoid  any  reflections  that  might  originate  from 
the  mesh  boundaries  and  affect  the  solution  later  in  time. 
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Figure  D5.  RAGE  3-D  problem  setup  for  GSA  configuration. 


Cross-section 

Driver  and  CVC  of  Striker  and  Flange 


The  walls  of  the  BLS  were  modeled  as  rigid  (non-responding)  materials. 
The  outer  shell  thickness  of  the  CYC  sections  was  artificially  increased  to  a 
thickness  of  0.79  in.,  about  1.6  times  the  true  thickness.  The  plot  in 
Figure  D5  suggests  that  the  shell  thickness  for  CVC-i  is  less  than  the 
others.  This  is  not  the  case;  rather  the  appearance  of  a  thinner  shell  is  a 
plotting  artifact  only.  The  shell  thickness  for  the  other  downstream 
sections  was  increased  even  more.  As  discussed  earlier,  it  was  necessary  to 
artificially  increase  the  shell  thickness  in  order  to  reduce  the  resolution 
requirements  and  ensure  tractable  run  times.  The  artificially  increased 
thickness  for  the  downstream  sections  is  not  an  issue,  as  they  contain  no 
vents;  however,  the  increased  thickness  of  the  CVC  shell  effectively 
reduces  the  vent  area  by  about  7.5  percent  on  average.  This  may  have  some 
effect  on  the  blow-down  later  in  time  and  hence  the  final  peak  impulse. 

The  surrounding  air  was  modeled  as  a  hydrodynamic  material.  The  air  was 
modeled  as  an  ideal  gas,  having  a  constant  gamma  of  1.4.  The  pressure  of 
the  driver  gas  was  1,298  psig.  It  was  assumed  to  be  at  an  ambient  temper¬ 
ature  of  77°F.  The  remainder  of  the  air  both  inside  and  outside  of  the  BLS 
was  set  to  ambient  conditions  (14.7  psia,  77°F).  More  detail  was  added  to 
the  mock-up  of  the  striker  mechanism,  with  the  legs  of  the  flange  included 
in  the  model.  Again,  no  attempt  was  made  to  include  the  grill  as  it  could 
not  be  adequately  resolved. 
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A  sensitivity  study  looking  at  refinement  level  was  initially  conducted  to 
assess  potential  bleed-through  in  the  CVC  section  for  the  3-D  problem 
setup.  The  findings  are  illustrated  in  Figure  D6,  which  depict  the  resultant 
velocity  field  at  2  msec  as  a  function  of  the  MLOR.  The  minimum  cell  size 
is  denoted  in  parentheses.  The  initial  cell  size  is  64  cm.  The  MLORs  con¬ 
sidered  were  6,  7,  and  8;  resulting  in  minimum  cell  sizes  of  2.0, 1.0,  and 
0.5  cm,  respectively.  It  is  clear  that  the  minimum  acceptable  resolution  for 
the  3-D  analysis  is  1.0  cm. 


Figure  D6.  Example  of  bleed-through  at  2  msec 
in  3-D  RAGE  analyses. 


(a)  MLOR  =  6  (2  cm) 


(b)  MLOR  =  7  (1  cm)  (c)  MLOR  =  8  (0.5  cm) 


Spatial  plots  of  the  resultant  velocity  field  at  select  times  are  provided  in 
Figure  D7.  The  scale  for  the  plots  on  the  right-hand  side  is  cast  in  terms  of 
cm/s.  The  scale  range  could  not  be  modified  using  the  post-processing 
tools  supplied  with  RAGE.  The  times  selected  for  the  plots  represent 
important  events  in  the  shock  propagation  problem.  By  5  msec,  the  shock 
front  has  traversed  the  CVC  section  and  has  entered  into  the  GSA  cone.  By 
10  msec,  the  front  has  entered  into  the  constant  area  section.  A  nearly 
planar  front  is  achieved  by  15  msec.  The  front  impinges  on  the  target  wall 
at  about  20  msec,  after  which  it  is  reflected.  The  reflected  shock  interacts 
with  any  residual  high  density  gas  propagating  downstream,  resulting  in  a 
highly  complex  shock  environment. 

Comparisons  between  the  measured  and  calculated  data  are  provided  in 
Figure  D8.  The  calculated  results  have  been  shifted  in  time  by  -3.15  msec 
to  match  the  measured  arrival  time.  The  time  shifting  was  done  to  ease 
comparisons  of  the  waveforms.  The  measured  data  were  fairly  consistent 
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(a)  5  msec 


(b)  10  msec 


(c)  15  msec 


(d)  20  msec 


(e)  25  msec 


(f)  30  msec 
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Figure  D8.  Time-shifted  gauge  comparisons  for  RAGE, 
GSA  configuration. 


(a)  Gauge  P3  (b)  Gauge  P4 


(c)  Gauge  PE 


across  the  target  wall,  indicating  that  the  incoming  shock  front  was  fairly 
planar  with  edge  effects  having  little  influence  on  the  reflected  pulse  (at 
least  over  the  time  frame  of  interest).  In  all  cases,  RAGE  significantly  over¬ 
predicts  the  peak  pressure  by  about  35  to  40  percent.  In  turn,  this  leads  to 
a  significant  overprediction  of  the  impulse.  The  degree  of  overprediction  is 
exaggerated  by  the  relatively  long  duration  of  the  initial  pulse.  As  noted  in 
the  2-DC  calculations,  there  is  a  significant  amount  of  ringing  in  the 
calculation.  The  ringing  diminishes  with  increasing  radial  offset  of  a 
gauge. 

It  is  worthwhile  to  compare  the  2-DC  and  3-D  RAGE  results.  Comparisons 
at  gauge  PE  are  provided  in  Figure  D9.  These  are  typical  of  those  at  the 
other  gauge  stations.  The  MLOR  was  held  fixed  for  the  duration  of  both 
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Figure  D9.  Comparisons  of  2-DC  and  3-D 
RAGE  caicuiations  for  gauge  PE. 
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the  2-DC  and  3-D  calculations.  Overall,  one  notes  excellent  agreement 
between  the  2-DC  and  3-D  calculations. 

The  3-D  rage  calculation  was  conducted  over  a  series  of  restarts  out  to  a 
final  analysis  time  of  52.6  msec.  The  output  from  each  calculation  was 
reviewed  to  obtain  timing  data.  The  calculation  was  run  at  a  MLOR  of  8 
(0.5  cm)  for  the  entire  calculation.  Over  the  duration  of  the  calculation,  the 
number  of  cells  in  the  problem  grew  from  13.4  to  279  million.  The  number 
of  processors  required  grew  from  128  at  startup  to  1,024  by  the  end  of  the 
calculation.  The  total  number  of  wall-clock  hours  to  conduct  the  entire 
calculation  was  429  hr.  It  took  approximately  282k  processor-CPU  hours 
for  the  calculation  to  reach  52.6  msec.  This  is  defined  as  the  CPU  time 
multiplied  by  the  number  of  processors,  with  the  quantity  summed  for  all 
restarts.  This  is  a  very  expensive  calculation. 

Calibrations  tests  with  8-ftx8-ft  C2SQ  configuration 

The  problem  setup  for  the  second  case  is  depicted  in  Figure  Dio.  This 
model  is  an  extension  of  the  GSA  setup  that  includes  the  C2SQ  and  SQi 
sections.  The  longitudinal  axis  of  the  BLS  ran  along  the  z-axis  in  the 
model,  with  the  origin  at  the  base  of  the  driver.  Half  symmetry  was 
assumed,  with  a  reflecting  plane  about  the  y-axis.  The  initial  cell  size  was 
64  cm.  The  number  of  level  zero  cells  along  the  x,  y,  and  z-axes  were  60, 
30,  and  40,  respectively.  Thus,  the  grid  encompassed  a  3,840-cm  x  1,920- 
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Figure  DIO.  RAGE  3-D  problem  setup  for  8-ftx8-ft 
C2SQ  configuration. 


cm  X  2560-cm  region  of  physical  space.  The  model  depicted  in  Figure  Dio 
was  built  up  with  overlaid  primitive  objects.  The  C2SQ  transition  section 
posed  some  challenges,  and  it  is  likely  that  the  curvature  in  the  corners  is 
sharper  than  indicated  in  the  photographs  of  Section  2.  A  CAD  model 
import  would  have  been  preferred;  however,  there  was  no  documentation 
(nor  example)  for  the  model  import. 

The  same  general  approach  for  the  problem  setup  from  the  first  case  was 
taken  here,  where  the  BLS  structure  was  treated  as  a  rigid  material,  and 
the  ambient  air  and  compressed  gas  were  modeled  as  an  ideal  gas  having  a 
ratio  of  specific  heats  of  1.4.  The  conditions  for  the  ambient  air  were  based 
on  an  average  of  that  specified  in  the  testing.  The  initial  driver  gas  pres¬ 
sure  was  set  to  1,298  psig  and  was  assumed  to  be  at  an  ambient  tempera¬ 
ture  of  6o°F.  The  remainder  of  the  air  both  inside  and  outside  of  the  BLS 
was  set  to  ambient  conditions  (14.7  psia,  6o°F).  As  with  the  previous  case, 
the  CVC  shell  thickness  has  been  artificially  increased  from  0.5  to  0.79  in. 
The  plot  in  Figure  Dii  suggests  that  the  shell  thickness  for  CVC-2  is  larger 
than  the  others.  This  is  not  the  case;  rather  the  appearance  of  a  thicker 
shell  is  a  plotting  artifact  only.  The  shell  thickness  for  the  other  down¬ 
stream  sections  was  increased  more.  This  is  most  obvious  in  the  C2SQ 
section.  This  transition  section  was  built  up  via  a  series  of  overlaid  primi¬ 
tive  objects.  The  thickness  was  exaggerated  to  avoid  gaps  being  induced  in 
the  model.  The  artificially  increased  thickness  of  the  downstream  sections 
should  have  no  effect  on  the  solution,  as  they  contain  no  vents;  however. 
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the  increased  thickness  of  the  CVC  shell  effectively  reduces  the  vent  area 
by  about  7.5  percent  on  average.  This  should  not  affect  the  initial  peak 
pressure,  but  could  potentially  affect  the  later  time  impulse  due  to  a 
reduction  in  the  rate  of  gas  venting. 


Spatial  plots  of  the  resultant  velocity  field  at  select  times  are  provided  in 
Figure  Dll.  The  scale  for  the  plots  on  the  right-hand  side  is  cast  in  terms  of 
cm/s.  The  scale  range  could  not  be  modified  using  the  post-processing 
tools  supplied  with  RAGE.  The  times  selected  for  the  plots  represent 
important  events  in  the  shock  propagation  problem.  By  5  msec,  the  shock 
front  has  traversed  the  CVC  section  and  has  entered  into  the  GSA  cone.  By 
10  msec,  the  front  has  entered  into  the  constant  area  section.  The  front 
impinges  on  the  target  wall  at  about  26  msec,  after  which  it  is  reflected. 

The  reflected  shock  interacts  with  any  residual  high-density  gas  propagat¬ 
ing  downstream,  resulting  in  a  highly  complex  shock  environment. 


Figure  Dll.  Velocity  field  from  RAGE  calculation, 
8-ftx8-ft  C2SQ  configuration. 


(a)  5  msec  (b)  10  msec 


(c)  25  msec 


(d)  35  msec 
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Comparisons  between  the  measured  and  calculated  data  are  provided  in 
Figure  D12  for  select  gauge  locations.  The  figure  contains  comparisons  for 
three  gauges  on  the  target  wall  (P6,  P5,  and  Pi)  along  with  the  side-on 
pressure  gauge  PR.  The  radial  distance  of  the  target  wall  gauges  from  the 
BLS  centerline  was  10.625,  24.5,  and  46.72  in.  for  gauges  P6,  P5,  and  Pi, 
respectively.  Thus,  the  gauges  selected  provide  a  good  indication  of  the 
influence  of  radial  offset  on  the  data.  It  is  important  to  note  that  the 
measured  data  were  fairly  consistent  across  the  target  wall,  indicating  that 
the  incoming  shock  front  was  fairly  planar  (at  least  over  the  time  frame  of 
interest).  The  calculated  pressure-time-histories  for  the  target  wall  gauges 
were  shifted  in  time  by  -2.93  msec  to  match  the  measured  arrival  time. 


Figure  D12.  Time-shifted  gauge  comparisons  for  RAGE, 
8-ftx8-ft  C2SQ  configuration. 


(a)  Gauge  PI 


(b)  Gauge  P5 


(d)  Gauge  PR 


(c)  Gauge  P6 
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Similarly,  the  calculated  side-on  gauge  data  were  shifted  by  -3.79  msec  to 
match  the  measured  arrival.  The  time  shifting  was  done  to  ease 
comparisons  of  the  waveform. 

RAGE  significantly  over-predicts  the  initial  peak  pressure  for  the  target 
wall  gauges  (on  the  order  of  40  to  45  percent).  In  turn,  this  leads  to  a 
significant  overprediction  in  the  impulse.  The  initial  side-on  pressure  at 
gauge  PR  is  also  overpredicted,  but  to  a  much  lesser  extent.  The  decay  rate 
behind  the  initial  peak  predicted  by  RAGE  is  higher  than  measured.  Fur¬ 
ther,  one  notes  a  longer  negative  phase  duration  in  the  RAGE  calculation 
as  compared  with  the  measured  data.  On  a  positive  note,  the  timing  of  the 
subsequent  spikes  is  captured  well  by  RAGE. 

The  RAGE  calculation  was  conducted  over  a  series  of  restarts  out  to  a  final 
analysis  time  of  57.2  msec.  The  output  from  each  calculation  was  reviewed 
to  obtain  timing  data.  The  calculation  was  run  at  a  MLOR  of  8  (0.5  cm)  for 
the  entire  calculation.  Over  the  duration  of  the  calculation,  the  number  of 
cells  in  the  problem  grew  from  38.5  to  509  million.  The  number  of  proces¬ 
sors  required  grew  from  256  at  startup  to  1,536  by  the  end  of  the  calcula¬ 
tion.  The  total  number  of  wall-clock  hours  to  carry  out  the  entire 
calculation  was  653  hr.  It  took  approximately  618k  processor-CPU  hours 
for  the  calculation  to  reach  57.2  msec.  This  is  defined  as  the  CPU  time 
multiplied  by  the  number  of  processors,  with  the  quantity  summed  for  all 
restarts.  This  is  also  a  very  expensive  calculation. 

BLS  8-ftx8-ft  C2SQ  configuration  with  box 

The  problem  setup  for  the  third  case  is  depicted  in  Figure  D13.  This  model 
is  identical  to  that  used  for  the  second  case,  with  the  exception  of  an 
instrumented  box  placed  in  the  C2SQ  section.  The  BLS  structure  and  box 
were  treated  as  rigid  materials;  the  ambient  air  and  compressed  gas  were 
modeled  as  an  ideal  gas  having  a  ratio  of  specific  heats  of  1.4.  The  condi¬ 
tions  for  the  ambient  air  were  based  on  an  average  of  that  specified  in  the 
testing.  The  initial  driver  gas  pressure  was  set  to  1,298  psig  and  was 
assumed  to  be  at  an  ambient  temperature  of  6o°F.  The  remainder  of  the 
air  both  inside  and  outside  of  the  BLS  was  set  to  ambient  conditions 
(14.7  psia,  6o°F).  As  with  the  previous  cases,  the  CVC  shell  thickness  was 
artificially  increased  from  0.5  to  0.79  in.  The  shell  thickness  for  the  other 
downstream  sections  was  increased  more.  This  is  most  obvious  in  the 
C2SQ  section.  This  transition  section  was  built  up  via  a  series  of  overlaid 
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Figure  D13.  RAGE  problem  setup  for  8-ftx8-ft 
C2SQ  configuration  with  box. 


Cross-Section 


primitive  objects.  A  box  was  inserted  into  the  model  along  the  base  of  the 
C2SQ  section.  Since  a  plane  of  symmetry  was  employed  in  the  problem 
setup,  it  was  only  possible  to  collect  information  on  one  side  of  the  box 
(right  face)  for  comparison  with  the  measured  data. 

Spatial  plots  of  the  resultant  velocity  field  at  selected  times  are  provided  in 
Figure  D14.  The  scale  for  the  plots  on  the  right-hand  side  is  cast  in  terms 
of  cm/s.  The  scale  range  could  not  be  modified  using  the  post-processing 
tools  supplied  with  RAGE.  Comparable  results  were  noted  with  the  previ¬ 
ous  case  for  times  leading  up  to  initial  interaction  with  the  box.  For  brev¬ 
ity,  they  have  not  been  included  here.  The  plots  in  Figure  D14  focus  on 
interaction  of  the  incident  and  reflected  shocks  with  the  box.  They  provide 
an  overall  view  of  the  BLS.  The  plots  in  Figures  D15  through  D17  provide  a 
close-up  view  of  the  region  about  the  box.  Independent  state  plots  of  the 
density,  pressure,  and  velocity  state  are  provided  at  select  times.  The  scale 
for  the  density  plots  is  in  g/cm3;  the  scale  for  the  pressure  plots  is  in 
dyne/cm2;  and  the  scale  for  the  velocity  magnitude  plots  is  in  cm/s.  It  is 
instructive  to  cross-reference  between  the  overall  and  close-up  views  of  the 
flow  state  when  analyzing  the  data. 
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Figure  D14.  RAGE  velocity  field,  8-ftx8-ft 
C2SQ  configuration  with  box. 


(e)  25  msec  (f)  28  msec 
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Figure  D14.  (Continued). 


(i)  32  msec  Q)  35  msec 


In  Figure  D14  at  20  msec,  the  incident  shock  front  is  about  to  impinge  on 
the  box.  By  22  msec,  the  front  has  traversed  the  box  and  is  diffracting 
about  the  back  side.  The  process  of  front  face  reflection,  passage  across  the 
top  surface,  and  diffraction  about  the  rear  of  the  box  is  best  illustrated  in 
the  detailed  state  plots.  After  traversing  the  box,  the  front  recombines 
downstream  to  form  an  almost  planar  front  again  by  25  msec.  It  reflects 
from  the  target  wall  at  about  26  msec.  The  reflected  shock  then  travels 
upstream  and  interacts  with  the  back  face  of  the  box  at  about  31  msec.  The 
reflected  shock  is  no  longer  planar  as  it  travels  toward  the  back  face  of  the 
box,  as  it  is  affected  by  the  residual  pressure  and  temperature  state  of  the 
air  in  the  region  between  the  box  and  target  wall.  By  35  msec,  the  reflected 
shock  has  traversed  the  box  and  passed  into  the  Cascade  section.  Although 
it  has  recovered  much  of  its  planar  structure,  this  will  be  short-lived  as  it 
continues  moving  upstream.  It  will  eventually  interact  with  the  residual 
gas  flow,  which  will  cause  it  to  lose  much  of  its  structure. 
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Figure  D15.  Flow  state  around  box  at  21  msec. 


(c)  Resultant  velocity 


Figure  D16.  Flow  state  around  box  at  22  msec. 


(a)  Density  (b)  Absolute  pressure 
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Figure  D16.  (Continued). 


(c)  Resuitant  veiocity 


Figure  D17.  Fiow  state  around  box  at  23  msec. 


(c)  Resuitant  veiocity 


Comparisons  between  the  measured  and  calculated  data  are  provided  in 
Figures  Di8  through  D21  for  the  gauges  positioned  on  the  box.  Compari¬ 
sons  for  the  gauges  located  on  the  target  wall  are  provided  in  Figure  D22. 
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Overall,  RAGE  tends  to  overpredict  the  peak  pressure  for  both  the  incident 
and  reflected  shocks,  which  in  turn,  leads  to  a  substantial  overprediction 
in  the  impulse.  This  is  readily  apparent  for  the  front  face  gauges  (see 
Figure  Di8),  with  less  of  a  disparity  noted  for  gauges  on  the  other  faces. 


Figure  D18.  Time-shifted  gauge  comparisons  for  RAGE,  front  face. 


(a)  Gauge  PBFl  (b)  Gauge  PBF3 


(c)  Gauge  PBF5 


(d)  Gauge  PBF6 


Impulse  (psi-msec) 


Overpressure  (psi)  Overpressure  (psi)  Overpressure  (psi) 
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Figure  D18.  (Continued). 


(e)  Gauge  PBF7 


(f)  Gauge  PBF8 


Figure  D19.  Time-shifted  gauge  comparisons  for  RAGE,  back  face. 


Time  (msec) 


Time  (msec) 


(a)  Gauge  PBBl 


(b)  Gauge  PBB2 


(c)  Gauge  PBB3 


(d)  Gauge  PBB4 


Overpressure  (psi)  Overpressure  (psi) 
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Figure  D20.  Time-shifted  gauge  comparisons  for  RAGE,  top  face. 


(a)  Gauge  PBTl 


(b)  Gauge  PBT2 


Time  (msec) 


(c)  Gauge  PBT3 


(d)  Gauge  PBT4 


(e)  Gauge  PBT5 


Impulse  (psi-msec)  Impulse  (psi-msec) 
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Figure  D21.  Time-shifted  gauge  comparisons  for  RAGE,  ieft  face. 


Time  (msec) 

(a)  Gauge  PBLl 


Time  (msec) 

(c)  Gauge  PBL3 


Time  (msec) 

(b)  Gauge  PBL2 


Time  (msec) 

(d)  Gauge  PBL4 


Figure  D22.  Time-shifted  gauge  comparisons  for  RAGE,  target  waii. 


Time  (msec) 


Time  (msec) 


(a)  Gauge  P3 


(b)  Gauge  P9 


Impulse  (psi-msec)  Impulse  (psi-msec) 
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Not  only  is  the  peak  pressure  overpredicted,  but  the  initial  rate  of  decay 
from  the  peak  is  generally  less  than  observed.  In  turn,  this  adds  to  the 
overprediction  noted  in  the  impulse,  as  this  is  simply  the  integration  of  the 
pressure-time-history.  Although  there  are  disparities  in  the  peaks  and  rate 
of  decay,  one  notes  reasonable  correlation  with  the  timing  of  the  peaks  and 
general  structure  of  the  waveforms.  For  example,  when  reviewing  the 
comparisons  for  the  front  face  gauges,  one  notes  two  distinct  peaks  in  the 
calculations  with  the  difference  in  arrival  of  the  two  peaks  correlating  well 
with  the  measured  data.  Similar  observations  apply  for  the  comparisons 
for  the  two  gauges  on  the  target  wall. 

The  RAGE  calculation  was  conducted  over  a  series  of  restarts  out  to  a  final 
analysis  time  of  6o  msec.  The  output  from  each  calculation  was  reviewed 
to  obtain  timing  data.  The  calculation  was  run  at  a  MLOR  of  8  (0.5  cm)  for 
the  entire  calculation.  Over  the  duration  of  the  calculation,  the  number  of 
cells  in  the  problem  grew  from  38.6  to  478  million.  The  number  of  proces¬ 
sors  required  grew  from  256  at  startup  to  2,048  by  the  end  of  the  calcu¬ 
lation.  The  total  number  of  wall-clock  hours  to  conduct  the  entire 
calculation  was  643  hr.  It  took  approximately  662k  processor-CPU  hours 
for  the  calculation  to  reach  60  msec.  This  is  defined  as  the  CPU  time  multi¬ 
plied  by  the  number  of  processors,  with  the  quantity  summed  for  all 
restarts.  This,  too,  is  a  very  expensive  calculation. 

Conclusions  and  recommendations 

Overall,  RAGE  tended  to  significantly  overpredict  the  initial  peak  pressure 
of  the  arriving  shock.  Furthermore,  the  initial  rate  of  decay  from  the  peak 
was  less  than  that  observed.  The  combined  effect  led  to  integrated 
impulses  that  were  significantly  greater  than  those  observed  (ranging  from 
20  to  45  percent  overprediction  depending  on  the  gauge).  Despite  the  dis¬ 
parities  in  the  peaks  and  rate  of  decay,  one  did  note  reasonable  correlation 
with  the  timing  of  the  peaks  and  general  structure  of  the  waveforms.  This 
was  most  evident  in  the  comparison  for  the  third  case,  where  multiple 
peaks  were  discerned  on  gauges  on  the  different  sides  of  the  box. 

The  reason  for  the  tendency  toward  overprediction  is  unclear.  The  early  2- 
DC  mesh  resolution  and  sensitivity  study  did  not  expose  any  “first  order” 
influences  on  the  calculations.  This  is  troubling  in  that  the  calculations 
appeared  relatively  insensitive  to  wide  variations  in  the  parameters  con¬ 
sidered.  It  is  conjectured  that  the  roots  of  the  problem  are  related  to  the 
use  of  an  ideal  gas  EOS.  The  SAGE  code  has  been  used  extensively  at 
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ERDC  for  modeling  airblast  resulting  from  explosive  detonations  (both 
ideal  and  non-ideal)  with  great  success.  The  initial  BLS  calculations  were 
conducted  using  SAGE  while  RAGE  was  being  installed  on  the  ERDC  HPC 
platform  Garnet.  The  results  for  the  two  codes  were  quite  similar,  with 
SAGE  also  exhibiting  a  tendency  for  overprediction  for  the  BLS  modeling. 
The  setup  for  the  early  SAGE  calculations  was  comparable  to  that  used  in 
the  RAGE  modeling.  It  also  treated  the  air  using  an  ideal  gas  EOS  with  a 
constant  gamma  of  1.4.  The  version  of  SAGE  in  use  at  ERDC  is  much 
older,  dating  back  to  the  2004  time  frame.  This  suggests  that  this  issue  of 
overprediction  may  be  unique  to  BLS-type  of  problems  involving  air  states 
in  relatively  low  pressure  regimes.  It  also  suggests  that  this  has  been  an 
unrecognized  issue  up  to  this  point. 

Of  course,  one  can  never  discount  the  potential  for  user  error.  Our  experi¬ 
ence  at  ERDC  has  been  exclusively  based  on  modeling  explosive 
detonations.  The  BLS  experiments  represent  a  vastly  different  class  of 
problems.  There  were  discussions  with  Leidos  regarding  the  BLS  problem 
setup.  Leidos  provided  many  constructive  suggestions,  in  particular,  focus¬ 
ing  on  alternative  options  in  the  hydrodynamic  treatment  as  well  as  alter¬ 
nate  setups  to  improve  run  times.  These  suggestions  were  generally  imple¬ 
mented  into  the  RAGE  modeling,  with  the  exception  of  modifying  the  BLS 
layout  to  ensure  edges  were  aligned  on  cell  boundaries.  This  really  is  not 
possible  given  that  many  BLS  sections  are  conical  in  nature.  Furthermore, 
one  cannot  align  the  circular  sections  in  a  3-D  setup.  For  these  reasons,  it 
was  decided  to  artificially  increase  the  BLS  shell  thickness  to  overcome  any 
difficulties  with  capturing  it  in  the  modeling. 

It  is  recommended  that  a  follow-on  analysis  be  conducted  with  the  LANL 
version  of  RAGE.  This  would  be  helpful  to  diagnose  the  issue  and  possibly 
pin-point  the  source  of  the  problem.  The  follow-on  work  could  be  very 
limited  in  scope.  A  simple  2-DC  calculation  would  suffice  and  would  pro¬ 
vide  a  quick  assessment  of  whether  the  issue  persists  in  the  LANL  version. 

It  is  also  recommended  that  the  user  manual  be  updated  to  reflect  the  cur¬ 
rent  version  of  the  code.  This  work  drew  upon  the  2005  user  manual  for 
SAGE,  which  is  only  an  incomplete  draft.  It  would  be  helpful  to  update  the 
user  documentation  and  develop  a  corresponding  set  of  example  prob¬ 
lems.  The  example  problems  would  illustrate  the  code’s  features  and  pro¬ 
vide  a  means  for  independently  testing  installations  on  new  HPC 
platforms. 
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Appendix  E;  SHAMRC  Modeling 

SHAMRC  Background 

The  Second-order  Hydrodynamic  Automatic  Mesh  Refinement  Code,  or 
SHAMRC,  is  a  government-owned  fluid  dynamics  suite  of  codes  developed 
and  maintained  by  Applied  Research  Associates  (ARA).  It  is  able  to  solve 
2-D  as  well  as  3-D  problems  using  a  finite-difference  scheme  and  can  be 
run  either  in  serial  or  parallel  processing  environments.  Some  of  the 
applications  for  which  the  code  is  used  include  high  explosive  detonations, 
nuclear  explosive  detonations,  structure  loading,  airblast,  thermal  effects, 
munitions  blast  and  fragmentation,  shock  tube  phenomenology,  dust  and 
debris  dispersion,  and  atmospheric  shock  propagation.  The  code  employs 
a  structured  Eulerian  grid  that  can  make  use  of  adaptive  mesh  refinement 
(AMR)  for  dividing  the  computational  domain  into  smaller  grids  at  several 
levels  of  refinement  to  provide  high-resolution  results.  Other  capabilities 
of  the  code  include  multiple  geometries,  non-responding  and  responding 
structures,  two-phase  flow  with  non-interactive  and  interactive  drag- 
sensitive  particles,  several  atmosphere  models,  multiple  materials,  a  large 
material  library,  HE  detonations,  a  K-e  turbulence  model,  water  and  dust 
vaporization,  and  a  predictive  metal  and  gaseous  afterburn  model. 
SHAMRC  is  second-order  accurate  in  space  and  time  and  is  fully  conser¬ 
vative  of  mass,  momentum,  and  energy. 

SHAMRC  is  better  described  as  a  suite  of  codes  as  it  is  composed  of  multi¬ 
ple  programs  that  are  used  for  setting  up,  running,  and  post-processing 
each  calculation  (Crepeau  et  al.  2012;  Happ  et  al.  2014).  The  main  prog¬ 
rams  that  are  used  for  achieving  this  are  listed  and  explained  here. 
Whenever  the  users  want  to  accomplish  any  of  these,  they  must  first  call 
SHAMRC  to  generate  the  source  code  (referred  to  as  “make”)  and  eventu¬ 
ally  execute  it  (simply  called  “run”).  This  makes  sure  that  the  code  is  tail¬ 
ored  and  optimized  specifically  for  the  problem  at  hand. 

1.  KEEL  -  generate  initial  computational  grid 

2.  SHARC  -  solver  that  advances  grid  conditions  in  time 

3.  PULL  -  plot  conditions  in  computational  grid  at  various  times 

4.  STASRT/STAPLT  -  process,  plot,  and  extract  station  (gauge)  data 

This  appendix  covers  all  of  the  calculations  that  were  performed  for  the 
three  different  BLS  configurations  using  SHAMRC. 
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Problem  setup 

All  of  the  SHAMRC  calculations  were  run  on  the  ERDC  HPC  platform 
Garnet.  Garnet  is  a  Cray  XE6  that  has  4,716  compute  nodes  with  32  cores 
per  node.  Each  compute  node  has  a  core  speed  of  2.5  GHz  and  64  GB  of 
accessible  memory.  Three  different  BLS  configurations  were  considered  in 
this  modeling  effort,  i.e.,  the  GSA  configuration,  the  8-ftx8-ft  C2SQ 
configuration,  and  the  8-ftx8-ft  C2SQ  configuration  with  box.  For 
simplicity,  these  will  be  referred  as  the  first,  second,  and  third  cases, 
respectively. 

Although  the  GSA  configuration  could  have  been  modeled  in  2-D  due  to  its 
axis-symmetric  geometry,  all  three  cases  were  modeled  in  3-D.  However, 
since  all  three  cases  were  symmetric  with  respect  to  the  BLS  centerline,  at 
least  a  one  half-symmetry  model  was  considered  for  each  case.  For  these 
half-symmetry  models,  reflective  boundary  conditions  were  specified  on 
the  plane  of  symmetry.  All  other  boundaries,  including  the  ones  for  the 
fully-3-D  models,  were  specified  as  transmissive,  i.e.,  flow  can  exit  the 
mesh  through  these  boundaries  but  it  is  not  allowed  to  flow  back  into  the 
grid. 

The  mesh  used  for  the  three  cases  was  flat,  i.e.,  it  did  not  involve  any  type 
of  AMR  nor  was  any  grading  applied.  Therefore,  the  cell  size  was  fixed 
throughout  the  computational  domain.  The  cell  size  was  selected  so  that  it 
would  be  less  than  the  smallest  dimension  in  the  problem.  Not  taking  into 
account  the  diaphragm-catch  grill  at  the  end  of  the  CVC-i  component,  the 
smallest  dimension  in  the  BLS  was  0.5  in.  (1.27  cm),  which  corresponded 
to  the  diameter  of  the  striker  as  well  as  the  thicknesses  of  the  outer  BLS 
shell  for  several  sections.  For  this  reason,  the  largest  cell  size  considered  in 
any  of  the  analyses  presented  herein  was  1.0  cm.  Other  resolutions  consid¬ 
ered  include  a  0.5-cm  cell  size  for  each  of  the  three  cases,  and  a  0.25-cm 
cell  size  for  the  case  having  the  box  in  the  flow. 

Table  Ei  summarizes  relevant  information  about  the  grid  resolution  and 
problem  size  for  each  case.  Additionally,  it  includes  a  “Nomenclature”  row 
used  for  identifying  each  model  considered.  Three  mesh  resolutions  were 
considered  and  are  referred  to  as  the  coarse,  medium,  and  fine  meshes. 
The  corresponding  cell  sizes  were  1.0,  0.5,  and  0.25  cm,  respectively. 
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Table  El.  Summary  of  models  considered  for  the 
three  BLS  configurations. 


Case 

1 

2 

3 

Nomenclature 

Coarse 

Medium 

(HS) 

Coarse 

Medium 

(HS) 

Coarse 

Medium 

Medium 

(HS) 

Fine 

(HS) 

Cell  size  (cm) 

1.0 

0.5 

1.0 

0.5 

1.0 

0.5 

0.5 

0.25 

Half-symmetry? 

No 

Yes 

No 

Yes 

No 

No 

Yes 

Yes 

Total  cells  (millions) 

256 

1024 

272 

1088 

272 

2176 

1088 

8704 

Models  that  make  use  of  half-symmetry  are  denoted  with  “(HS).”  Half¬ 
symmetry  was  not  used  in  any  of  the  coarse-meshed  models,  as 
computational  time  was  manageable,  and  results  could  be  obtained  in  a 
relatively  short  amount  of  time.  For  problems  having  smaller  cell  sizes, 
half-symmetry  was  essential  for  reducing  the  computational  cost. 

Only  one  material  was  specified  in  the  input  files  of  all  models,  i.e.,  air.  All 
of  the  solid  components  of  the  BLS,  including  the  walls,  were  defined  as 
“islands.”  An  island,  in  SHAMRC  jargon,  refer  to  a  non- responding 
structure  in  the  mesh.  These  rigid  structures  have  no  momentum  or 
energy  associated  with  them,  and  their  density  is  equal  to  one.  Because  of 
these  properties,  all  of  the  solid  components  that  were  modeled  were 
entirely  reflective.  Additionally,  all  island/solid  geometries  were 
constructed  by  using  3-D  simple  entities,  such  as  cylinders,  rectangular 
cuboids,  cones,  and  tetrahedrons. 

For  air,  two  different  initial  conditions  were  specified  in  all  cases.  First,  a 
constant  density  and  pressure  atmosphere  model  was  set  for  the  ambient 
air,  which  used  a  density  of  1.2250-3  g/cm3  and  an  internal  energy  of 
2.04409  ergs/g.  Secondly,  a  pressure  of  8.9493907  dynes/cm^  (1,298  psi) 
and  a  density  of  1.0820-1  g/cm3  were  set  for  the  compressed  air  in  the 
driver.  The  density  corresponded  to  air  at  a  pressure  of  1,298  psi  and  a 
temperature  of  66°F.  Additionally,  three  different  equations  of  state  were 
considered,  a  Doan-Nickel  version,  a  fixed  gamma  law  (ideal  gas),  and  a 
tabular  air  from  Los  Alamos  National  Laboratory  (SESAME  tables).  A 
comparison  between  these  is  made  in  this  appendix  for  the  BLS  having  a 
GSA-configuration  and  using  a  cell  resolution  of  1.0  cm.  Based  on 
recommendations  by  an  experienced  SHAMRC  user,  the  Doan-Nickel  EOS 
was  selected  for  use  on  all  cases. 
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Tracers,  or  stations  in  SHAMRC  jargon,  were  placed  throughout  the  mesh 
at  the  specified  gauge  locations  and  then  adjusted  slightly  to  account  for 
mesh  resolution.  This  was  done  to  avoid  issues  with  mixed  cells  being 
partially  filled  with  air  and  “islands.”  As  a  general  rule,  stations  residing  on 
a  wall  were  offset  by  1.5  times  the  resolution  level  of  each  model.  For 
meshes  composed  of  1.0,  0.5,  and  0.25-cm  cells,  the  spatial  offsets  were 
1-55  0.75,  and  0.375  cm,  respectively.  This  condition  positioned  stations  at 
least  one  cell  away  from  their  corresponding  wall.  All  tracers  were  fixed  in 
space,  with  material  allowed  to  impinge  upon  and  flow  through  their 
location. 

Time  offset 

When  compared  with  the  measured  data,  all  of  the  incident  pressure 
pulses  in  the  calculated  results  show  a  lag  in  the  time  of  arrival,  usually 
about  3.0  to  4.0  msec.  The  reason  for  the  difference  in  arrival  times  is 
unknown  at  this  time.  To  ease  comparisons,  the  calculated  arrival  times 
are  shifted  throughout  to  match  the  measured  data.  The  amount  of  time 
offset  is  also  affected  by  the  mesh  resolution.  The  time  offsets  applied  for 
each  calculation  are  noted  in  the  ensuing  discussion. 

EOS  selection 

SHAMRC  has  several  EOS  options.  Of  these,  only  three  appeared  appro¬ 
priate  for  modeling  the  air  in  the  BLS:  (a)  fixed  gamma  law,  (b)  Doan- 
Nickel  real  air,  and  (c)  tabular  air  from  LANE  SESAME  tables.  The  first 
one,  the  fixed  gamma  law,  is  commonly  known  as  the  ideal  gas  law.  This 
EOS  assumes  a  constant  gamma  of  1.4  for  air.  The  second  one,  the  Doan- 
Nickel,  was  developed  by  Larry  Doan  and  George  Nickel  (1963)  and  is 
based  on  the  fact  the  gamma  value  of  air,  which  is  the  ratio  between  the 
specific  heat  at  constant  pressure  and  specific  heat  at  constant  volume,  will 
vary  from  a  value  of  1.4  when  incident  pressures  go  beyond  300  psi.  The 
third  one  is  based  on  the  SESAME  EOS  library  developed  by  Los  Alamos 
National  Laboratory. 

These  three  EOSs  were  tested  using  the  BLS  GSA  configuration,  which  is 
the  first  case  out  of  the  three  different  experimental  configurations.  The 
mesh  size  for  these  numerical  tests  was  fixed  at  1.0  cm.  Figure  Ei  com¬ 
pares  the  results  obtained  from  these  EOSs  for  three  gauges  located  on  the 
BLS  target  wall:  P3,  P4,  and  PE.  As  seen  in  the  figure,  there  is  no  signifi¬ 
cant  difference  in  the  shape  of  the  calculated  pressure  pulse,  although 
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Figure  El.  Influence  of  EOS  on  calculated  result,  GSA  configuration. 


(a)  Gauge  P3  (b)  Gauge  P4 


(c)  Gauge  PE 

there  are  slight  differences  in  the  magnitude  of  the  pressure  at  specific 
time-steps.  These  slight  differences  in  the  pressure  are  enough  to  make  the 
integrated  impulses  deviate  from  each  other  after  about  6o  msec,  which  is 
the  time  at  which  the  pressure  wave  reflects  from  the  target  wall  for  the 
second  time.  However,  the  deviations  do  not  exceed  3  percent  by  the  time 
the  numerical  analysis  was  concluded  at  100  msec.  It  is  interesting  to  note 
that  there  is  no  consistency  in  regard  to  which  EOS  will  produce  the  high¬ 
est  impulse  at  the  end  of  the  pressure-time  records.  For  gauge  P3,  the 
highest  impulse  is  attributed  to  the  Doan-Nickel  EOS,  while  the  lowest  one 
corresponds  to  the  ideal  gas  EOS.  This  trend  is  exactly  the  opposite  for 
gauge  PE.  For  gauge  P4,  all  three  EOS  produce  almost  the  same  integrated 
impulse.  Another  interesting  observation  is  presented  in  Figure  E2,  which 
shows  zoomed  views  of  the  initial  peak  pressures  calculated  for  Gauges  P3, 
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Figure  E2.  Calculated  initial  pressure  peaks  as  a  function  of  EOS. 


(a)  Gauge  P3  (b)  Gauge  P4 


P4,  and  PE.  Note  that  in  all  three  records,  the  EOS  having  the  earliest  time 
of  arrival  is  the  ideal  gas  followed  by  the  SESAME  and  Doan-Nickel  EOS 
(in  that  order).  Moreover,  a  similar  behavior  can  be  seen  for  the  initial 
pressure  peak.  In  this  case,  the  highest  pressure  is  always  produced  by  the 
SESAME  EOS,  followed  by  the  ideal  gas  and  then  the  Doan-Nickel  EOS. 

Comparisons  with  the  measured  data  are  provided  in  Figure  E3.  At  first,  it 
is  obvious  that,  although  the  shape  of  the  pressure-time  records  agree,  the 
initial  calculated  pressure  peak  and  impulse  are  both  overpredicted  for 
each  EOS.  Following  the  initial  pressure  pulse,  during  the  pressure  decay 
and  second  reflection,  the  computed  pressures  agree  with  the  measured 
data.  Since  none  of  the  EOS  provided  a  best  estimate  to  the  measured 
data,  additional  criteria  were  considered  to  select  one  EOS  to  be 
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Figure  E3.  EOS  comparisons  with  measured  data,  GSA  configuration. 


Time  (msec)  Time  (msec) 

(a)  Gauge  P3  (b)  Gauge  P4 


Time  (msec) 


(c)  Gauge  PE 

used  in  the  validation  calculations.  Considering  that  the  driver  pressure  of 
1,298  psi  will  produce  an  incident  pressure  above  300  psi  close  to  the 
driver,  the  Doan-Nickel  EOS  seemed  like  a  good  choice.  Moreover,  the  use 
of  this  EOS  was  recommended  by  an  experienced  SHAMRC  user. 

GSA  configuration,  caiibration  Tests  30a  and  30b 

The  first  of  the  three  BLS  configurations  that  was  tested  was  the  GSA 
configuration.  Two  numerical  models  were  run  in  SHAMRC  for  this  case, 
with  comparisons  made  with  measurements  on  the  target  wall.  The  time 
offsets  applied  to  the  calculated  data  for  each  model  considered  are 
outlined  in  Table  E2.  Time  offsets  were  applied  to  each  gauge.  The  time 
offsets  for  the  two  models  are  comparable. 
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Table  E2.  Time  offsets  for  GSA  configuration  models  (in  msec). 


Gauge 

Coarse 

Medium  (HS) 

Gauge 

Coarse 

Medium  (HS) 

PI 

-3.64 

-3.64 

P6 

-3.83 

-3.80 

P2 

-3.83 

-3.77 

P7 

-3.86 

-3.85 

P3 

-3.77 

-3.76 

P9 

-3.79 

-3.79 

P4 

-3.83 

-3.79 

PIO 

-3.82 

-3.82 

P5 

-3.79 

-3.79 

PE 

-3.86 

-3.84 

With  regard  to  the  computational  resources  that  were  used  for  both  mod¬ 
els,  the  coarse  model  was  run  for  31.0  hr  using  64  processors  while  the 
half-symmetric  “medium”  model  was  started  on  128  processors  and  even¬ 
tually  increased  to  256  processors.  The  128  processors  were  used  for 
23.0  hr  while  the  256  processors  were  used  for  70.5  hr.  Both  models  were 
run  up  to  no  msec. 

Comparisons  between  the  measured  and  calculated  data  are  provided  in 
Figure  E4.  Disregarding  the  magnitudes  of  the  calculated  pressure  and 
integrated  impulse,  the  overall  shape  of  the  pressure-time  records  agree. 
The  general  shape  of  the  waveform  is  composed  of  an  initial  peak  that 
occurs  at  about  16  to  17  msec,  and  a  second  peak  of  less  amplitude  that 
occurs  around  60  msec.  These  two  pressure  spikes  cause  the  integrated 
impulse  to  be  composed  of  two  noticeable  “bumps.” 

In  general,  the  initial  calculated  pressure  peak  is  higher  than  that  mea¬ 
sured  by  as  much  as  50  percent.  This  overprediction  causes  the  impulse 
from  the  numerical  results  to  be  greater  than  that  measured.  During  the 
first  pressure  decay,  at  usually  about  30  to  40  msec,  the  calculated 
pressure  briefly  goes  below  the  experimental  values  and  then  returns  to 
similar  values.  This  underprediction  is  seen  as  well  during  the  pressure 
decay  following  the  second  reflected  pressure  pulse.  Another  observation 
seen  during  the  second  pressure  pulse  of  some  records  (e.g.,  P3,  P4,  and 
P6)  is  that  multiple  pressure  reflections  in  one  of  the  calculated  results 
cause  the  integrated  impulse  to  deviate  from  the  other  calculated  one. 
There  is  no  consistency  in  regards  to  which  of  the  two  models  will  end  up 
having  the  highest  impulse.  However,  there  is  a  tendency  that  the  impulse 
of  the  numerical  results  will  be  overpredicted  by  as  much  as  10  to 
20  percent.  An  exception  to  this  trend  is  gauge  P2  (Figure  E4(b)).  The 
initial  calculated  peak  pressure  correlates  well  with  that  measured, 
resulting  in  excellent  comparisons  with  the  integrated  impulse. 
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Figure  E4.  Time-shifted  gauge  comparisons  for 
SHAMRC,  GSA  configuration. 


(a)  Gauge  PI 


(b)  Gauge  P2 


Time  (msec) 


Time  (msec) 


(c)  Gauge  P3 


(d)  Gauge  P4 


Time  (msec) 


Time  (msec) 


(e)  Gauge  P5 


(f)  Gauge  P6 
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Figure  E4.  (Continued). 


(g)  Gauge  P7  (h)  Gauge  P9 


(i)  Gauge  PIO  Q)  Gauge  PE 

A  series  of  pressure  contours  are  presented  in  Figure  E5  showing  the  pres¬ 
sure  wave  as  it  traverses  the  BLS  and  reflects  from  the  target  wall.  The 
legend  at  the  top  of  the  figure  shows  the  pressure  magnitude  in  dynes/cm^. 
Note  that  the  pressure  front  becomes  planar  once  it  leaves  the  conical  sec¬ 
tion.  However,  a  mismatch  in  the  radial  dimensions  of  the  Cascade  sec¬ 
tions  seen  at  20  msec  makes  the  portion  of  the  pressure  front  farther  from 
the  BLS  centerline  lag.  This,  in  turn,  affects  the  measured  TOA  for  gauges 
having  the  largest  radial  distance  from  the  BLS  centerline. 
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Figure  E5.  Pressure  contours  from  SHARMC  calculation, 
GSA  configuration. 
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7.239E+06  9.653E+06 

dynes/sq-cm 


(a)  0.0  msec 


(b)  7.93  msec 


(c)  10.0  msec 


(d)  14.56  msec 


(e)  20.0  msec 
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Figure  E5.  (Continued). 


(g)  40.0  msec 


Calibration  tests  with  8-ftx8-ft  C2SQ  configuration 

Similar  to  the  GSA  configuration,  two  models  were  run  for  this  second 
case.  The  first  was  a  fully-3-D  mesh  with  a  1.0-cm  resolution.  The  second 
involved  a  half-symmetry  mesh  with  a  0.5-cm  resolution.  In  the  tables  and 
figures  of  this  section,  these  two  runs  are  identified  as  “Coarse”  and 
“Medium  (HS)”,  respectively.  For  the  coarse  model,  256  processors  were 
used  for  a  total  wall  time  of  about  26.6  hr.  For  the  half-symmetry  medium 
model,  256  processors  were  used  for  about  68.9  hr. 

As  mentioned  before,  the  numerical  results  obtained  for  this  case  had  a 
later  time  of  arrival  than  that  measured.  Because  of  this,  the  calculated 
pressure-time-histories  were  shifted  to  match  the  measured  values. 

Table  E3  outlines  the  time  shifts  applied  to  each  gauge.  In  general,  the 
offsets  for  the  coarse  models  are  less  than  those  for  the  medium  one.  The 
reason  for  this,  as  mentioned  in  the  setup  section  of  this  appendix,  is  that 
the  stations  in  the  coarser  model  were  spatially  offset  twice  the  distance 
(1.5  vs.  0.75  cm)  than  those  in  the  medium  model.  This  positioned  the 
gauges  in  the  coarse  model  closer  to  the  driver  resulting  in  slightly  earlier 
times  of  arrival. 
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Table  E3.  Time  offsets  for  8  ftx8  ft  C2SQ  configuration 
models  (in  msec). 


Gauge 

Coarse 

Medium  (HS) 

Gauge 

Coarse 

Medium  (HS) 

PI 

-3.43 

-3.95 

P8 

-3.36 

-3.87 

P2 

-3.30 

-3.83 

P9 

-3.40 

-3.93 

P3 

-3.32 

-3.81 

PIO 

-3.57 

-4.07 

P4 

-3.39 

-3.91 

Pll 

-3.57 

-4.08 

P5 

-3.36 

-3.88 

P12 

-3.43 

-3.96 

P6 

-3.35 

-3.87 

PR 

-3.87 

-4.06 

P7 

-3.37 

-3.87 

Comparisons  between  the  calculated  and  measured  data  are  provided  in 
Figure  E6.  Except  for  a  few  gauges,  these  comparisons  are  not  very  differ¬ 
ent  from  the  results  that  were  presented  for  the  first  case.  First,  the  shape 
of  the  pressure-time  records  is  dominated  by  two  pressure  pulses,  which 
correspond  to  the  first  and  second  interaction  of  the  pressure  wave  with 
the  gauges  on  the  target  wall.  Second,  the  two  pressure  pulses  are  “echoed” 
in  the  integrated  impulses  as  the  shape  of  these  can  be  explained  as  having 
two  curvilinear  sections.  Third,  the  initial  peak  pressure  tends  to  be  over¬ 
predicted,  causing  the  impulse  to  be  overpredicted  as  well.  A  noticeable 
difference  with  respect  to  the  first  case  is  that  the  calculated  pressures  fol¬ 
lowing  the  first  pressure  decay  are  in  better  agreement  with  the  experi¬ 
mental  data. 

There  were  two  gauges  on  the  target  wall  whose  behavior  differed  mark¬ 
edly  from  the  other  gauges.  The  measured  data  for  gauge  P7  is  suspiciously 
low  as  compared  to  the  other  gauges,  suggesting  an  anomalous  measure¬ 
ment  at  this  gauge  location.  For  gauge  P8,  the  calculated  pressure  and 
impulse  matches  the  measured  data  until  the  first  pressure  decay,  where 
the  pressure  is  slightly  lower,  resulting  in  an  underprediction  of  the 
impulse. 
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Figure  E6.  Time-shifted  gauge  comparisons  for  SHAMRC, 
8-ftx8-ft  C2SQ  configuration. 
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Figure  E6.  (Continued). 


(g)  Gauge  P7 


(h)  Gauge  P8 


(i)  Gauge  P9 


0)  Gauge  PIO 


(k)  Gauge  Pll 


(i)  Gauge  P12 
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Comparisons  between  the  calculated  and  measured  data  for  the  side-on 
pressure  gauge  PE  are  provided  in  Figure  E7.  The  same  trends  noted  pre¬ 
viously  apply  here.  That  is,  there  is  an  overall  tendency  of  overprediction 
in  the  calculated  results.  Given  that  this  is  a  side-on  gauge,  there  is  a 
noticeable  difference  in  the  shape  and  frequency  of  secondary  reflections 
when  compared  to  the  other  gauges.  The  calculation  appears  to  capture 
the  timing  of  these  secondary  reflections;  however,  there  are  notable 
differences  in  the  amplitudes  of  the  peaks. 


Figure  E7.  Comparisons  for  side-on 
pressure  gauge  PE. 


Figure  E8  shows  pressure  contours  of  the  compressed  air  front  traveling 
from  the  driver  to  the  target  wall  of  the  BLS.  The  contours  at  time  zero  are 
not  included,  as  they  were  not  radically  different  from  what  was  shown  for 
the  first  case.  One  main  difference  between  these  contours  and  the  ones 
shown  earlier  is  that  it  takes  more  time  for  the  pressure  front  to  reach  the 
end-wall,  given  that  the  overall  length  of  the  BLS  is  greater  for  this  case. 
Another  difference  is  the  inclusion  of  the  grill  at  the  location  of  the  first 
vent.  However,  it  is  not  clear  from  these  images  if  the  grill  induces  signifi¬ 
cant  turbulence  into  the  flow  as  it  traverses  the  grill. 
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Figure  E8.  Pressure  contours  from  SHAMRC  calculation, 
8-ftx8-ft  C2SQ  configuration. 
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BLS  8-ftx8-ft  C2SQ  configuration  with  box 

The  third  case  had  the  same  configuration  as  the  second  with  the  only 
difference  being  an  instrumented  box  placed  in  the  C2SQ  section.  Four 
models  were  considered  for  this  case.  They  included  two  fully-3-D  models 
having  resolutions  of  1.0  and  0.5  cm,  and  two  half-symmetry  models  hav¬ 
ing  resolutions  of  0.5  and  0.25  cm.  Each  of  these  models  required  varying 
computational  resources.  The  first  model  (coarse,  fully-3-D)  utilized 
128  processors  for  a  total  of  52  hr  of  wall  clock  time.  The  second  model 
(medium,  fully-3-D)  was  initially  run  for  23.5  hr  using  128  processors,  but 
the  number  of  processors  had  to  be  increased  to  decrease  the  computa¬ 
tional  time.  The  model  was  restarted  using  256  processors  for  47  hr,  after 
which  it  was  restarted  again  using  512  processors  for  48  hr.  The  third 
model  (medium,  half-symmetry)  was  started  on  128  processors  for 
47.75  hr  and  restarted  with  256  processors  for  53.5  hr.  The  fourth  model 
(fine,  half-symmetry)  was  started  on  4,096  processors  for  176.3  hr,  which 
advanced  the  simulation  to  about  45.0  msec.  The  number  of  processors 
was  doubled  to  8,192  for  about  42.1  hr,  which  progressed  the  simulation 
time  to  58.58  msec.  Finally,  the  number  of  processors  was  increased  to 
16,384  for  29.5  hr,  which  carried  the  calculation  out  to  a  little  over 
60.0  msec.  Given  the  excessive  cost  of  the  fine  meshed  calculation,  it  was 
decided  not  to  continue  it  any  further. 

The  times  for  which  the  numerical  results  were  shifted  to  match  the  mea¬ 
sured  TOA  are  summarized  in  Tables  E4  and  E5.  There  is  no  consistency 
in  which  of  the  four  models  required  the  largest  time  offset,  but  the  one 
that  required  the  smallest  offsets  was  the  coarse  model.  Moreover,  the 
medium  fully-3-D  and  half-symmetric  models  required  identical  time  off¬ 
sets  for  most  of  gauges,  which  was  expected  given  that  the  cell  sizes  of  both 
problems  were  equal.  Gauges  PBRi  and  PBR2  do  not  have  a  time  offset 
recorded  for  the  medium  model,  nor  are  they  included  in  the  pressure¬ 
time  plots  because  they  were  incorrectly  positioned  in  the  model. 
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Table  E4.  Time  offsets  for  gauges  on  the  target  wall  (in  msec). 


Gauge 

Coarse 

Medium 

Medium  (HS) 

Fine  (HS) 

P3 

-2.97 

-3.22 

-3.20 

-3.14 

P9 

-3.06 

-3.30 

-3.30 

-3.23 

Table  E5.  Time  offsets  for  gauges  on  box  (in  msec). 


Gauge 

Coarse 

Medium 

Medium  (HS) 

Fine  (HS) 

PBBl 

-2.97 

-3.47 

-3.48 

-3.42 

PBB2 

-3.06 

-3.48 

-3.50 

-3.44 

PBB3 

-3.24 

-3.49 

-3.51 

-3.44 

PBB4 

-3.31 

-3.51 

-3.51 

-3.45 

PBFl 

-3.30 

-3.49 

-3.49 

-3.44 

PBF3 

-3.32 

-3.50 

-3.50 

-3.47 

PBF5 

-3.28 

-3.50 

-3.50 

-3.47 

PBF6 

-3.31 

-3.49 

-3.49 

-3.45 

PBF7 

-3.39 

-3.49 

-3.52 

-3.48 

PBF8 

-3.30 

-3.50 

-3.50 

-3.47 

PBLl 

-3.30 

-3.51 

-3.51 

-3.47 

PBL2 

-3.32 

-3.50 

-3.48 

3.44 

PBL3 

-3.30 

-3.50 

-3.50 

-3.44 

PBL4 

-3.29 

-3.49 

-3.49 

-3.42 

PBRl 

-3.29 

N/A 

-3.48 

-3.42 

PBR2 

-3.30 

N/A 

-3.50 

-3.45 

PBR3 

-3.30 

-3.48 

-3.48 

-3.48 

PBR4 

-3.31 

-3.50 

-3.50 

-3.45 

PBTl 
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-3.52 
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Comparisons  between  the  calculated  and  measured  data  are  provided  in 
Figures  E9  through  E14.  The  first  set  of  records  show  the  gauges  located 
on  the  target  wall  (see  Figure  E9).  The  behavior  seen  in  these  plots  is  not 
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too  different  from  that  of  the  cases  discussed  earlier,  where  the  initial  peak 
pressure  is  overpredicted  in  the  calculations,  causing  the  integrated 
impulse  to  be  higher  than  that  measured. 

The  comparisons  for  gauges  located  on  the  back  face  of  the  box  are  pro¬ 
vided  in  Figure  Eio.  The  shapes  of  these  records  differ  from  the  ones 
shown  for  gages  P3  and  P9  as  they  are  facing  downstream,  away  from  the 
driver.  The  first  peak  observed  in  the  records  is  attributed  to  the  first  inter¬ 
action  of  the  pressure  wave  as  it  travels  over  the  box.  The  second  is  attrib¬ 
uted  to  the  reflected  wave  that  is  travelling  upstream.  Overall,  the 
pressures  and  impulses  are  overpredicted  in  the  calculation.  An  exception 
to  this  is  gauge  PBB3  in  which  the  duration  of  the  first  pressure  pulse  in 
the  calculation  is  shorter  than  that  measured,  resulting  in  better  correla¬ 
tion  with  the  measured  impulse.  Additionally,  gauge  PBB4  displays  the 
biggest  difference  in  the  integrated  impulse  due  to  a  high  peak  pressure 
exhibited  in  the  second  pressure  pulse  in  the  plot. 


The  comparisons  for  gauges  located  on  the  front  face  of  the  box  are  pro¬ 
vided  in  Figure  Ell.  These  results  are  reminiscent  of  all  the  gauges  previ¬ 
ously  shown  that  were  positioned  on  the  target  wall.  That  is,  almost  all 
these  records  show  overpredictions  in  the  calculated  pressures  as  well  as 
in  the  integrated  impulse.  The  obvious  exception  to  this  overprediction  is 
gauge  PBE6,  whose  measured  impulse  is  11  to  12  percent  higher  than  the 
calculated  results.  This  is  caused  by  a  lower  calculated  pressure  associated 
with  the  second,  reflected  pressure  pulse. 


(a)  Gauge  P3 


(b)  Gauge  P9 


Impulse  (psi-msec) 


Overpressure  (psi)  Overpressure  (psi) 
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Figure  ElO.  Time-shifted  gauge  comparisons  for  SHAMRC,  back  face. 
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Figure  Ell.  Time-shifted  gauge  comparisons  for  SHAMRC,  front  face. 
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The  comparisons  for  gauges  on  the  left  and  right  side  of  the  box  are  shown 
in  Figures  E12  and  E13,  respectively.  Both  sets  of  results  are  not  very  dif¬ 
ferent  from  each  other.  They  all  contain  an  initial  pressure  pulse  of  short 
duration  followed  by  a  second  and  larger  pressure  pulse  that  lasts  longer. 
A  third,  less  intense  pressure  pulse  can  be  seen  at  about  80  msec. 


Figure  E12.  Time-shifted  gauge  comparisons  for  SHAMRC,  ieftface. 
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(c)  Gauge  PBL3 


(d)  Gauge  PBL4 
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Figure  E13.  Time-shifted  gauge  comparisons  for  SH AM RC,  right  face. 
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Overall,  the  calculated  results  are  in  good  agreement  with  the  experimen¬ 
tal  data  for  gauges  on  the  right  and  left  side  of  the  box,  with  the  greatest 
disparities  noted  for  the  coarse-meshed  model.  On  most  pressure-time- 
histories  shown  in  Figures  Ei2  and  E13,  the  final  impulse  for  the  coarse 
model  is  underpredicted  by  5  to  10  percent.  This  is  mostly  the  result  of 
underpredicting  the  peak  pressure  during  the  second  pressure  pulse. 
Recall  that  gauges  in  the  models  are  offset  1.5  cells  away  from  the  BLS 
component  in  which  they  are  located  to  avoid  issues  with  mixed  cells,  and 
that  the  model  having  the  largest  offset  is  the  one  with  the  largest  cell  size 
(i.e.,  the  coarse  model).  It  is  likely  that  the  gauges  for  the  coarser  model 
residing  on  both  sides  of  the  box  experienced  less  reflected  pressure  than 
the  gauges  in  finer  models,  resulting  in  the  data  for  these  gauges  being  less 
than  that  measured. 
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The  comparisons  for  gauges  located  on  the  top  face  of  the  box  are  provided 
in  Figure  E14.  These  plots  have  a  shape  that  is  similar  to  the  ones  shown 
for  gauges  on  the  left  and  right  side  of  the  box.  Moreover,  the  calculated 
results  compare  well  with  the  measured  data.  The  peak  pressures  are 
almost  identical,  and  the  duration  of  each  pressure  pulse  is  comparable. 
The  measured  impulse  is  bounded  by  the  calculated  ones  on  every  plot.  On 
three  out  of  the  five  plots,  the  coarse  model  exhibits  the  lowest  impulse.  It 
is  likely  that  this  is  related  to  the  spatial  offset  of  gauges  in  this  model  from 
the  box,  as  it  is  larger  for  the  coarse  model  than  for  any  other  one. 

Figure  E15  shows  a  top  and  side  view  of  pressure  contours  displaying  the 
compressed  air  wave  interacting  with  the  box.  There  are  three  relevant 
phenomena  that  can  be  observed  in  these  images.  First,  one  can  see  how 
the  pressure  front  becomes  almost  planar  as  it  nears  the  box.  This  planar¬ 
ity,  however,  does  not  translate  into  a  constant  magnitude  of  pressure  as 
several  shades  of  blue  can  be  clearly  seen  inside  of  the  compressed  air, 
indicating  different  values  of  pressure.  Secondly,  a  green  intense  color  can 
be  seen  on  the  front  of  the  structure  in  the  21.2-  to  22.0-msec  frames.  This 
is  caused  by  the  wave  reflecting  from  the  box,  resulting  in  an  increase  to 
the  pressure  magnitude.  Finally,  several  low  pressure  regions  can  be 
observed  on  both  sides  of  the  box  as  well  as  above  it  as  the  pressure  front 
engulfs  the  box.  These  dark  blue  spots  are  an  indication  of  vortices  occur¬ 
ring  as  the  pressure  wave  impinges  against  the  box  and  diffracts  around  it. 


Overpressure  (psi)  Overpressure  (psi) 
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Figure  E14.  Time-shifted  gauge  comparisons  for  SHAMRC,  top  face. 
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Figure  E15.  Pressure  contours  around  box  in  SHARMC  caicuiation. 
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Figure  E15.  (Continued). 
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Other  spatial  images 

Various  images  are  presented  in  this  section  showing  how  several  BLS 
components  were  actually  represented  in  models  having  different  mesh 
resolution.  The  first  image,  shown  in  Figure  E16,  shows  a  cross-section 
view  of  the  striker  flange  at  the  mid-section  of  the  flange.  The  curvature  of 
the  flange  structure  is  better  represented  as  the  mesh  is  refined.  SHAMRC 
has  an  option  termed  “shores,”  which  can  be  used  to  alleviate  these  “stair¬ 
step”  structures.  Shores  are  essentially  cells  composed  evenly  of  an  island 
and  a  fluid.  In  this  case,  a  shore  boundary  would  pass  through  the  diago¬ 
nal  of  a  cell  for  a  more  accurate  representation  of  the  structure.  However, 
they  are  only  available  in  2-D  calculations,  and  their  use  was  not 
employed. 
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Figure  E16.  Striker  flange  cross-section. 


As  seen  in  the  previous  figure,  the  striker  flange  is  composed  of  four  bars 
that  attach  the  outer  ring  to  a  smaller  disk  at  the  center  of  the  flange. 
Attached  to  the  front  and  back  of  this  disk  are  the  striker  and  spring  mech¬ 
anism,  respectively.  These  were  modeled  as  two  cylinders  of  different 
diameters.  Figure  E17  shows  a  top  view  at  half-width  of  the  striker  flange 
that  was  shown  in  the  previous  figure.  The  horizontal  sections  in  these 
images  are  the  two  horizontal  bars  shown  in  the  previous  figure.  The  verti¬ 
cal  section  above  the  bars  seen  in  the  three  images  is  the  spring  mecha¬ 
nism.  The  striker  is  below  the  bars  and  can  only  be  seen  in  the  last  two 
images.  As  with  the  previous  figure.  Figure  E17  shows  how  this  component 
was  actually  represented  in  models  having  different  resolutions.  In  the 
first  image,  one  can  see  that  the  striker  does  not  even  show  up  in  the 
model.  Moreover,  the  thickness  of  the  BLS  wall  around  the  striker  flange 
appears  to  be  greater  than  the  one  shown  for  the  other  two  models.  In  the 
center  image,  the  striker  is  visible  and  the  thickness  of  the  BLS  wall 
appears  to  be  less  than  in  the  first  image.  In  the  image  on  the  right,  one 
can  actually  notice  a  difference  in  the  thickness  and  radius  of  the  CYC-i 
component  inside  the  Ci  section  (above  the  horizontal  bars  of  the  striker 
flange)  with  respect  to  the  striker  flange  and  driver  section  (around  and 
below  the  horizontal  bars  of  the  striker  flange). 

The  grill,  as  it  was  represented  in  models  having  different  resolution,  is 
shown  in  Figure  E18.  The  grill  is  better  resolved  as  the  mesh  resolution 
increases.  For  the  model  having  a  cell  size  of  1.0  cm,  the  grill  is  poorly 
replicated.  When  the  cell  size  is  decreased  by  half,  some  of  the  features  of 
the  grill  start  showing  up.  However,  at  0.5  cm,  the  grill  is  still  missing  six 
of  the  horizontal  bars.  It  is  not  until  the  cells  size  is  decreased  to  0.25  cm 
that  the  grill  is  captured  in  its  entirety. 
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Figure  E17.  Top  view  of  striker  and  fiange. 
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Figure  E18.  Griii  cross-section. 


(a)  1.0  cm  (b)  0.5  cm  (c)  0.25  cm 


Figure  E19  shows  profiles  at  mid-width  of  the  BLS  for  each  of  the  cases 
considered.  Several  components  can  be  seen  in  the  images,  including  the 
striker,  the  striker  flange,  the  spring  mechanism,  the  grill  and  the 
instrumented  box. 

Issues 

Several  issues  were  noted  while  running  SHAMRC.  Two  of  these  issues  are 
related  to  generating  plots  using  PULL,  while  the  third  one  is  about  a 
warning  message  that  kept  showing  up  while  the  calculations  were 
ongoing. 
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Figure  E19.  Cross-section  of  the  BLS  for  the 
three  cases  considered. 


The  first  issue  worth  mentioning  is  about  generating  spatial  plots  at  spe¬ 
cific  times  of  the  simulation  using  PULL.  Like  many  other  CFD  codes,  spa¬ 
tial  plots  can  be  created  from  times  at  which  the  code  generated  a  “dump” 
file.  A  dump  file  contains  all  of  the  necessary  information  to  restart  the 
simulation.  The  file  size  for  these  dump  files  can  be  quite  large,  depending 
on  the  size  of  the  computational  grid.  For  the  calculations  shown  in  this 
appendix,  these  dump  files  were  around  a  dozen  gigabytes  for  the  coarsest 
model  to  several  hundred  gigabytes  for  the  finest  one.  So  if  one,  as  an 
example,  wanted  to  create  plots  at  lOO  specific  time-steps,  one  would  first 
need  to  generate  dump  files  at  each  of  the  time-steps.  This  will  be  a  prob¬ 
lem  if  one  has  limited  disk  space.  Assuming  that  the  size  of  a  dump  file  is 
around  lOO  gigabytes,  generating  these  lOO  spatial  plots  would  require 
10,000  gigabytes  of  space. 

The  second  issue  worth  mentioning  is  related  to  generating  spatial  plots 
for  large  problems  (problems  having  more  than  i  billion  cells).  In  this 
case,  the  time  needed  to  generate  one  spatial  image  is  proportional  to  the 
number  of  cells  in  the  problem.  This  is  not  an  issue  for  problems  having  a 
relatively  small  amount  of  cells,  but  it  is  definitely  perceptible  for  larger 
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problems.  As  an  example,  consider  the  half-symmetry  fine-meshed  model 
which  consists  of  about  8.7  billion  cells.  It  took  the  code  about  79  min  to 
generate  the  single  spatial  image  depicted  in  Figure  E19.  So  if  one  wanted 
to  generate  20  spatial  plots,  it  would  take  more  than  a  full  day  to  generate 
them.  It  is  relevant  to  mention  that  the  program  PULL,  used  to  generate 
spatial  images,  can  only  be  run  in  serial. 

The  third  and  final  issue  worth  mentioning  is  about  a  message  that  kept 
showing  up  in  the  output  files  of  the  models  considered.  The  message 
warned  about  the  energy  being  below  1.0E6  ergs  in  some  areas  of  the  mesh 
being  used  by  a  particular  processor.  Below  is  an  example  of  such 
message: 


On  processor  19: 

Energy  below  l.e6  in  MOVEY  at  i=  1  j=  162  k=  139 
h: 

O.OOOOOE+00 -8.61164E-02  7.23796E+04  2.23971E+09  1.53152E-04  1.01315E+06 
Energy  below  location  (2.500000E-01  -1.975000E+01  1.687500E+02). 


The  user’s  manual  warns  that  this  message  is  printed  whenever  the  inter¬ 
nal  energy  of  a  cell  drops  below  1.0E6  ergs  in  the  advection  phase.  The 
code  reacts  by  determining  the  energy  deficit  and  reducing  the  kinetic 
energy  of  the  cell  by  the  deficit  amount  so  that  the  total  energy  of  the  cell 
will  be  1.0E6.  According  to  the  manual,  this  generally  happens  when  the 
stability  factor  is  too  high  for  the  calculation.  Therefore,  it  was  decided  to 
choose  stability  factors  that  matched  what  was  recommended  by  the  user’s 
manual.  For  the  multiple  3-D  simulations  shown  in  this  appendix,  an 
initial  value  of  0.1  was  chosen  as  the  stability  factor  for  the  calculation  of 
the  Courant  condition  time-step.  This  factor  was  allowed  to  gradually 
increase  by  0.02  up  to  a  maximum  of  0.5. 
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